Bioinformatics systems, apparatuses, and methods executed on an integrated circuit processing platform

ABSTRACT

A system, method and apparatus for executing a sequence analysis pipeline on genetic sequence data includes an integrated circuit formed of a set of hardwired digital logic circuits that are interconnected by physical electrical interconnects. One of the physical electrical interconnects forms an input to the integrated circuit connected with an electronic data source for receiving reads of genomic data. The hardwired digital logic circuits are arranged as a set of processing engines, each processing engine being formed of a subset of the hardwired digital logic circuits to perform one or more steps in the sequence analysis pipeline on the reads of genomic data. Each subset of the hardwired digital logic circuits is formed in a wired configuration to perform the one or more steps in the sequence analysis pipeline.

CROSS REFERENCE TO RELATED APPLICATIONS

This application is a continuation of U.S. patent application Ser. No.14/988,666 entitled “Bioinformatics Systems, Apparatuses, and MethodsExecuted on an Integrated Circuit Processing Platform,” filed Jan. 5,2016; U.S. patent application Ser. No. 14/988,666 is a continuation ofU.S. patent application Ser. No. 14/284,307, entitled “BioinformaticsSystems, Apparatuses, and Methods Executed on an Integrated CircuitProcessing Platform,” filed May 21, 2015; U.S. patent application Ser.No. 14/284,307 is a continuation of U.S. patent application Ser. No.14/279,063, entitled “Bioinformatics Systems, Apparatuses, and MethodsExecuted on an Integrated Circuit Processing Platform,” filed May 15,2014, a continuation in part of: U.S. patent application Ser. No.14/180,248, entitled “Bioinformatics Systems, Apparatuses, and MethodsExecuted on an Integrated Circuit Processing Platform,” filed Feb. 13,2014, now patented as U.S. Pat. No. 9,014,989, and a continuation ofU.S. patent application Ser. No. 14/158,758, entitled “BioinformaticsSystems, Apparatuses, and Methods Executed on an Integrated CircuitProcessing Platform,” filed Jan. 17, 2014; U.S. patent application Ser.No. 14/279,063 is a continuation in part of U.S. patent application Ser.No. 14/180,248, now patented as U.S. Pat. No. 9,014,989, a continuationin part of U.S. patent application Ser. No. 14/179,513, entitled“Bioinformatics Systems, Apparatuses, and Methods Executed on anIntegrated Circuit Processing Platform,” filed Feb. 12, 2014, acontinuation in part of U.S. patent application Ser. No. 14/158,758, andclaims the benefit of and priority to under 35 U.S.C. 119(e) of U.S.Provisional Application Ser. No. 61/753,775, titled, “System and Methodfor Bioinformatics Processor,” filed Jan. 17, 2013; U.S. ProvisionalApplication Ser. No. 61/822,101, titled, “Bioinformatics ProcessorPipeline Based on Population Inference,” filed May 10, 2013, f U.S.Provisional Application Ser. No. 61/823,824, titled, “BioinformaticsProcessing System,” filed May 15, 2013; U.S. Provisional ApplicationSer. No. 61/826,381 titled, “System and Method for Computation GenomicsPipeline,” filed May 22, 2013, U.S. Provisional Application Ser. No.61/910,868, titled, “Bio-Informatics Systems and Methods Executed On aHardware Processing Platform,” filed Dec. 2, 2013, U.S. ProvisionalApplication Ser. No. 61/988,128 titled, “Bioinformatics Systems,Apparatuses, and Methods Executed on an Integrated Circuit ProcessingPlatform,” filed May 2, 2014, U.S. Provisional Application Ser. No.61/984,663 titled, “Bioinformatics Systems, Apparatuses, and MethodsExecuted on an Integrated Circuit Processing Platform” filed Apr. 25,2014, and U.S. Provisional Application Ser. No. 61/943,870 titled,“Dynamic Genome Reference Generation for Improved NGS Accuracy andReproducibility” filed Feb. 24, 2014; U.S. patent application Ser. No.14/179,513 is a continuation of U.S. patent application Ser. No.14/158,758; U.S. patent application Ser. No. 14/158,758 claims thebenefit of and priority under 35 U.S.C. 119(e) of U.S. ProvisionalApplication Ser. No. 61/753,775, U.S. Provisional Application Ser. No.61/822,101, U.S. Provisional Application Ser. No. 61/823,824, U.S.Provisional Application Ser. No. 61/826,381, U.S. ProvisionalApplication Ser. No. 61/910,868, U.S. Provisional Application Ser. No.61/988,128, U.S. Provisional Application Ser. No. 61/984,663, and, U.S.Provisional Application Ser. No. 61/943,870; U.S. patent applicationSer. No. 14/180,248, entitled “Bioinformatics Systems, Apparatuses, andMethods Executed on an Integrated Circuit Processing Platform,” filedFeb. 13, 2014, now patented as U.S. Pat. No. 9,014,989 is a continuationin part of Ser. No. 14/158,758, entitled “Bioinformatics Systems,Apparatuses, and Methods Executed on an Integrated Circuit ProcessingPlatform,” filed Jan. 17, 2014. The disclosures of the above-identifiedpatent applications are hereby incorporated by reference in theirentirety.

TECHNICAL FIELD

The subject matter described herein relates to bioinformatics, and moreparticularly to systems, apparatuses, and methods for implementingbioinformatic protocols, such as performing one or more functions foranalyzing genomic data on an integrated circuit, such as on a hardwareprocessing platform.

BACKGROUND

A goal for health care researchers and practitioners is to improve thesafety, quality, and effectiveness of health care for every patient.Personalized health care is directed to achieving these goals on anindividual level. For instance, “genomics” and/or “bioinformatics” arefields of study that aim to facilitate the safety, the quality, and theeffectiveness of prophylactic and therapeutic treatments on apersonalized, individual level. Accordingly, by employing genomicsand/or bioinformatics techniques, the identity of an individual'sgenetic makeup, e.g., his or hers genes, may be determined and thatknowledge may be used in the development of therapeutic and/orprophylactic regimens, including drug treatments, that are personalizedto the individual, thus, enabling medicine to be tailored to meet eachperson's individual needs.

The desire to provide personalized care to individuals is transformingthe health care system. This transformation of the health care system islikely to be powered by breakthrough innovations at the intersection ofmedical science and information technology such as is represented by thefields of genomics and bioinformatics. Accordingly, genomics andbioinformatics are key foundations upon which this future will be built.Science has evolved dramatically since the first human genome was fullysequenced in 2000 at a total cost of over $1 Billion. Today, we are onthe verge of high resolution sequencing at a cost of less than $1K pergenome, making it economically feasible for the first time to move outof the research lab and into widespread adoption for medical care.Genomic data, therefore, may become a vital input to diagnosticscreening, therapeutic and/or prophylactic drug discovery, and/ordisease treatment.

More particularly, genomics and bioinformatics are fields concerned withthe application of information technology and computer science to thefield of molecular biology. In particular, bioinformatics techniques canbe applied to process and analyze various genomic data, such as from anindividual so as to determine qualitative and quantitative informationabout that data that can then be used by various practitioners in thedevelopment of prophylactic and therapeutic methods for preventing or atleast ameliorating diseased states, and thus, improving the safety,quality, and effectiveness of health care on an individualized level.

Because of its focus on advancing personalized healthcare,bioinformatics, therefore, promotes individualized healthcare that isproactive, instead of reactive, and this gives the patient theopportunity to become more involved in their own wellness. Typically,this can be achieved through two guiding principles. First, federalleadership can be provided to support research that addresses theseindividual aspects of disease and disease prevention, such as with theultimate goal of shaping diagnostic and preventative care to match eachperson's unique genetic characteristics. Additionally, a “network ofnetworks” may be created to aggregate health care data to helpresearchers establish patterns and identify genetic “definitions” toexisting diseases.

An advantage of employing bioinformatics technologies in such instancesis that the qualitative and/or quantitative analyses of molecularbiological data can be performed on a broader range of sample sets at amuch higher rate of speed and often times more accurately, thusexpediting the emergence of a personalized healthcare system.

Accordingly, in various instances, the molecular data to be processed ina bioinformatics based platform typically concerns genomic data, such asDeoxyribonucleic acid (DNA) data. For example, a well-known method forgenerating DNA data involves DNA sequencing. DNA sequencing can beperformed manually, such as in a lab, or may be performed by anautomated sequencer, such as at a core sequencing facility, for thepurpose of determining the genetic makeup of a sample of an individual'sDNA. The person's genetic information may then be used in comparison toa referent so as to determine its variance therefrom. Such variantinformation may then be subjected to further processing and used todetermine or predict the occurrence of a diseased state in theindividual.

For instance, manual or automated DNA sequencing may be employed todetermine the sequence of nucleotide bases in a sample of DNA, such as asample obtained from a subject. Using various different bioinformaticstechniques these sequences may then be strung together to generate thegenomic sequence of the subject. This sequence may then be compared to areference genomic sequence to determine how the genomic sequence of thesubject varies from that of the reference. Such a process involvesdetermining the variants in the sampled sequence and presents a centralchallenge to bioinformatics methodologies.

For example, a central challenge in DNA sequencing is buildingfull-length genomic sequences, e.g., chromosomal sequences, from asample of genetic material that can be compared to a reference genomicsequence such as to determine the variants in the sampled full-lengthgenomic sequences. In particular, the methods employed in sequencingprotocols do not produce full-length chromosomal sequences of the sampleDNA.

Rather, sequence fragments, typically from 100-1,000 nucleotides inlength, are produced without any indication as to where in the genomethey align. Therefore, in order to generate full length chromosomalgenomic constructs, these fragments of DNA sequences need to be mapped,aligned, merged, and/or compared to a reference genomic sequence.Through such processes the variants of the sample genomic sequences fromthe reference genomic sequences may be determined.

However, as the human genome is comprised of approximately 3.1 billionbase pairs, and as each sequence fragment is typically only from 100 to500 nucleotides in length, the time and effort that goes into buildingsuch full length genomic sequences and determining the variants thereinis quite extensive often requiring the use of several different computerresources applying several different algorithms over prolonged periodsof time.

In a particular instance, thousands to millions of fragments of DNAsequences are generated, aligned, and merged in order to construct agenomic sequence that approximates a chromosome in length. A step inthis process may include comparing the DNA fragments to a referencesequence to determine where in the genome the fragments align.

A number of such steps are involved in building chromosome lengthsequences and in determining the variants of the sampled sequence.Accordingly, a wide variety of methods have been developed forperforming these steps. For instance, there exist commonly used softwareimplementations for performing one or a series of such steps in abioinformatics system. However, a common characteristic of such softwarebased bioinformatics methods and systems is that they are laborintensive, take a long time to execute on general purpose processors,and are prone to errors.

A bioinformatics system, therefore, that could perform the algorithmsimplemented by such software in a less labor and/or processing intensivemanner with a greater percentage accuracy would be useful. However, evenas we approach the “$1000 Genome”, the cost of analyzing, storing andsharing this raw digital data has far outpaced the cost of producing it.This data analysis bottleneck is a key obstacle standing between theseever-growing raw data and the real medical insight we seek from it.

Accordingly, presented herein are systems, apparatuses, and methods forimplementing a genomics and/or bioinformatic protocols, such as forperforming one or more functions for analyzing genomic data, forinstance, on an integrated circuit, such as on a hardware processingplatform. For example, as set forth herein below, in variousimplementations, a hardware accelerator, such as an integrated circuit,may be employed in performing such bioinformatics related tasks wherethe integrated circuit may be formed of one or more hardwired digitallogic circuits, which may be interconnected by a plurality of physicalelectrical interconnects, that can be arranged as a set of processingengines, wherein each processing engine is capable of being configuredto perform one or more steps in a bioinformatics genetic analysisprotocol. An advantage of this arrangement is that the bioinformaticsrelated tasks may be performed in a manner that is faster than thesoftware typically engaged for performing such tasks. Such hardwareaccelerator technology, however, is currently not typically employed inthe genomics and/or bioinformatics space.

SUMMARY

This present disclosure is related to performing a task such as in abioinformatics protocol. In various instances, a plurality of tasks areperformed, and in some instances these tasks are performed in a mannerso as to form a pipeline, wherein each task and/or its substantialcompletion acts as a building block for each subsequent task until adesired end result is achieved. Accordingly, in various embodiments, thepresent disclosure is directed to performing one or more methods on oneor more apparatuses wherein the apparatus has been optimized forperforming those methods. In certain embodiments, the one or moremethods and/or one or more apparatuses are formulated into one or moresystems.

For instance, in certain aspects, the present disclosure is directed tosystems, apparatuses, and methods for implementing genomics and/orbioinformatic protocols such as, in various instances, for performingone or more functions for analyzing genetic data on an integratedcircuit, such as implemented in a hardware processing platform. Forexample, in one aspect, a bioinformatics system is provided. The systemmay involve the performance of various bioanalytical functions that havebeen optimized so as to be performed faster and/or with increasedaccuracy. The methods for performing these functions may be implementedin software or hardware solutions. Accordingly, in certain instances,methods are presented where the method involves the performance of analgorithm where the algorithm has been optimized in accordance with themanner in which it is to be implemented. In particular, where thealgorithm is to be implemented in a software solution, the algorithmand/or its attendant processes, has been optimized so as to be performedfaster and/or with better accuracy for execution by that media.Likewise, where the functions of algorithm are to be implemented in ahardware solution, the hardware has been designed to perform thesefunctions and/or their attendant processes in an optimized manner so asto be performed faster and/or with better accuracy for execution by thatmedia.

Accordingly, in one aspect, presented herein are systems, apparatuses,and methods for implementing bioinformatic protocols, such as forperforming one or more functions for analyzing genetic data, forinstance, via one or more optimized algorithms and/or on one or moreoptimized integrated circuits, such as on one or more hardwareprocessing platforms. Hence, in one instance, methods are provided forimplementing one or more algorithms for the performance of one or moresteps for analyzing genomic data in a bioinformatics protocol. Inanother instance, methods are provided for implementing the functions ofone or more algorithms for the performance of one or more steps foranalyzing genomic data in a bioinformatics protocol, wherein thefunctions are implemented on an integrated circuit formed of one or morehardwired digital logic circuits. In such an instance, the hardwireddigital logic circuits may be interconnected, such as by one or aplurality of physical electrical interconnects, and may be arranged tofunction as one or more processing engines. In various instances, aplurality of hardwired digital logic circuits are provided, whichhardwired digital logic circuits are configured as a set of processingengines, wherein each processing engine is capable of performing one ormore steps in a bioinformatics genetic analysis protocol.

More particularly, in one instance, a system for executing a sequenceanalysis pipeline such as on genetic sequence data is provided. Thesystem may include one or more of an electronic data source, a memory,and an integrated circuit. For instance, in one embodiment, anelectronic data source is included, where in the electronic data sourcemay be configured for providing one or more digital signals, such as adigital signal representing one or more reads of genetic data, forexample, where each read of genomic data includes a sequence ofnucleotides. Further, the memory may be configured for storing one ormore genetic reference sequences, and may further be configured forstoring an index, such as an index of the one or more genetic referencesequences.

Further still, the integrated circuit may be formed of a set ofhardwired digital logic circuits such as where the hardwired digitallogic circuits are interconnected, e.g., by a plurality of physicalelectrical interconnects. In various instances, one or more of theplurality of physical electrical interconnects may include an input,such as to the integrated circuit, and may further be connected with theelectronic data source, so as to be able to receive the one or morereads of genomic data. In various embodiments, the hardwired digitallogic circuits may be arranged as a set of processing engines, such aswhere each processing engine is formed of a subset of the hardwireddigital logic circuits, and is configured so as to perform one or moresteps in the sequence analysis pipeline, such as on the plurality ofreads of genomic data. In such instances, each subset of the hardwireddigital logic circuits may be in a wired configuration so as to performthe one or more steps in the sequence analysis pipeline.

Accordingly, in various instances, a plurality of hardwired digitallogic circuits are provided wherein the hardwired digital logic circuitsare arranged as a set of processing engines, wherein one or more of theprocessing engines may include one or more of a mapping module and/or analignment module and/or a sorting module. For instance, in variousembodiments, the one or more of the processing engines may include amapping module, which mapping module may be in a wired configuration andfurther be configured for accessing the index of the one or more geneticreference sequences from the memory, such as by one or more of theplurality of physical electronic interconnects, for example, so as tomap the plurality of reads to one or more segments of the one or moregenetic reference sequences.

Additionally, in various embodiments, the one or more of the processingengines may include an alignment module, which alignment module may bein the wired configuration and may be configured for accessing the oneor more genetic reference sequences from the memory, such as by one ormore of the plurality of physical electronic interconnects, for example,so as to align the plurality of reads to the one or more segments of theone or more genetic reference sequences. Further, in variousembodiments, the one or more of the processing engines may include asorting module, which sorting module may be in the wired configurationand may be configured for accessing the one or more aligned reads fromthe memory, such as by one or more of the plurality of physicalelectronic interconnects, for example, so as to sort each aligned read,such as according to its one or more positions in the one or moregenetic reference sequences. In such instances, the one or more of theplurality of physical electrical interconnects may include an outputfrom the integrated circuit, such as for communicating result data fromthe mapping module and/or the alignment module and/or the sortingmodule.

In various instances, the integrated circuit may include a mastercontroller so as to establish the wired configuration for each subset ofthe hardwired digital logic circuits, for instance, for performing theone or more of mapping, aligning, and/or sorting, which functions may beconfigured as one or steps in a sequence analysis pipeline. Further, invarious embodiments, the integrated circuit may be configured as a fieldprogrammable gate array (FPGA) having hardwired digital logic circuits,such as where the wired configuration may be established uponmanufacture of the integrated circuit, and thus may be non-volatile. Inother various embodiments, the integrated circuit may be configured asan application specific integrated circuit (ASIC) having hardwireddigital logic circuits.

In certain instances, the integrated circuit and/or the memory may behoused on an expansion card, such as a peripheral component interconnect(PCI) card, for instance, in various embodiments, the integrated circuitmay be a chip having a PCIe card. In various instances, the integratedcircuit and/or chip may be a component within a sequencer, such as anautomated sequencer, and/or in other embodiments, the integrated circuitand/or expansion card may be accessible via the internet, e.g., cloud.Further, in some instances, the memory may be a volatile random accessmemory (RAM).

Accordingly, in one aspect, an apparatus for executing one or more stepsof a sequence analysis pipeline, such as on genetic data, is providedwherein the genetic data includes one or more of a genetic referencesequence(s), an index of the one or more genetic reference sequence(s),and/or a plurality of reads, such as of genetic data. In variousinstances, the apparatus may include an integrated circuit, whichintegrated circuit may include one or more, e.g., a set, of hardwireddigital logic circuits, wherein the set of hardwired digital logiccircuits may be interconnected, such as by one or a plurality ofphysical electrical interconnects. In certain instances, the one or moreof the plurality of physical electrical interconnects may include aninput, such as for receiving the plurality of reads of genomic data.Additionally, the set of hardwired digital logic circuits may further bein a wired configuration, so as to access the index of the one or moregenetic reference sequences, via one of the plurality of physicalelectrical interconnects, and to map the plurality of reads to one ormore segments of the one or more genetic reference sequences, such asaccording to the index.

In various embodiments, the index may include one or more hash tables,such as a primary and/or secondary hash table. For instance, a primaryhash table may be included, wherein in such an instance, the set ofhardwired digital logic circuits may be configured to do one or more of:extracting one or more seeds of genetic data from the plurality of readsof genetic data; executing a primary hash function, such as on the oneor more seeds of genetic data so as to generate a lookup address foreach of the one or more seeds; and accessing the primary hash tableusing the lookup address so as to provide a location in the one or moregenetic reference sequences for each of the one or more seeds of geneticdata. In various instances, the one or more seeds of genetic data mayhave a fixed number of nucleotides.

Further, in various embodiments, the index may include a secondary hashtable, such as where the set of hardwired digital logic circuits isconfigured for at least one of extending at least one of the one or moreseeds with additional neighboring nucleotides, so as to produce at leastone extended seed of genetic data; executing a hash function, e.g., asecondary hash function, on the at least one extended seed of geneticdata, so as to generate a second lookup address for the at least oneextended seed; and accessing the secondary hash table, e.g., using thesecond lookup address, so as to provide a location in the one or moregenetic reference sequences for each of the at least one extended seedof genetic data. In various instances, the secondary hash function maybe executed by the set of hardwired digital logic circuits, such as whenthe primary hash table returns an extend record instructing the set ofhardwired digital logic circuits to extend the at least one of the oneor more seeds with the additional neighboring nucleotides. In certaininstances, the extend record may specify the number of additionalneighboring nucleotides by which the at least one or more seeds isextended, and/or the manner in which the seed is to be extended, e.g.,equally by an even number of “x” nucleotides to each end of the seed.

Additionally, in one aspect, an apparatus for executing one or moresteps of a sequence analysis pipeline on genetic sequence data isprovided, wherein the genetic sequence data includes one or more of oneor a plurality of genetic reference sequences, an index of the one ormore genetic reference sequences, and a plurality of reads of genomicdata. In various instances, the apparatus may include an integratedcircuit, which integrated circuit may include one or more, e.g., a set,of hardwired digital logic circuits, wherein the set of hardwireddigital logic circuits may be interconnected, such as by one or aplurality of physical electrical interconnects. In certain instances,the one or more of the plurality of physical electrical interconnectsmay include an input, such as for receiving the plurality of reads ofgenomic data. Additionally, the set of hardwired digital logic circuitsmay further be in a wired configuration, so as to access the one or moregenetic reference sequences, via one of the plurality of physicalelectrical interconnects, to receive location information specifying oneor more segments of the one or more reference sequences, and to alignthe plurality of reads to the one or more segments of the one or moregenetic reference sequences.

In various instances, the wired configuration of the set of hardwireddigital logic circuits, are configured to align the plurality of readsto the one or more segments of the one or more genetic referencesequences, and further include a wave front processor that me be formedof the wired configuration of the set of hardwired digital logiccircuits. In certain embodiments, the wave front processor may beconfigured to process an array of cells of an alignment matrix, such asa matrix defined by a subset of the set of hardwired digital logiccircuits. For instance, in certain instances, the alignment matrix maydefine a first axis, e.g., representing one of the plurality of reads,and a second axis, e.g., representing one of the segments of the one ormore genetic reference sequences. In such an instance, the wave frontprocessor may be configured to generate a wave front pattern of cellsthat extend across the array of cells from the first axis to the secondaxis; and may further be configured to generate a score, such as foreach cell in the wave front pattern of cells, which score may representthe degree of matching of the one of the plurality of reads and the oneof the segments of the one or more genetic reference sequences.

In such an instance, the wave front processor may further be configuredso as to steer the wave front pattern of cells over the alignment matrixsuch that the highest score may be centered on the wave front pattern ofcells. Additionally, in various embodiments, the wave front processormay further be configured to backtrace one or more, e.g., all, thepositions in the scored wave front pattern of cells through previouspositions in the alignment matrix; track one or more, e.g., all, of thebacktraced paths until a convergence is generated; and generate a CIGARstring based on the backtrace from the convergence.

In certain embodiments, the wired configuration of the set of hardwireddigital logic circuits to align the plurality of reads to the one ormore segments of the one or more genetic reference sequences may includea wired configuration to implement a Smith-Waterman and/orBurrows-Wheeler scoring algorithm. In such an instance, theSmith-Waterman and/or Burrows-Wheeler scoring algorithm may beconfigured to implement a scoring parameter that is sensitive to basequality scores. Further, in certain embodiments, the Smith-Watermanscoring algorithm may be an affine Smith-Waterman scoring algorithm.

Accordingly, in one aspect, a method for executing a sequence analysispipeline such as on genetic sequence data is provided. The genetic datamay include one or more genetic reference sequences, one or more indexesof the one or more genetic reference sequences, and/or a plurality ofreads of genomic data. The method may include one or more of receiving,accessing, mapping, aligning, and/or sorting various iterations of thegenetic sequence data. For instance, in certain embodiments, the methodmay include receiving, on an input to an integrated circuit from anelectronic data source, one or more of a plurality of reads of genomicdata, wherein each read of genomic data may include a sequence ofnucleotides. In such an instance, the integrated circuit may be formedof a set of hardwired digital logic circuits such as are interconnectedby a plurality of physical electrical interconnects, which physicalelectrical interconnects may include one or more of the plurality ofphysical electrical interconnects comprising the input.

The method may further include accessing, by the integrated circuit onone or more of the plurality of physical electrical interconnects from amemory, the index of the one or more genetic reference sequences. Insuch an instance the method may include mapping, by a first subset ofthe hardwired digital logic circuits of the integrated circuit, theplurality of reads to one or more segments of the one or more geneticreference sequences. Additionally, the method may include accessing, bythe integrated circuit on one or more of the plurality of physicalelectrical interconnects from the memory, the one or more geneticreference sequences; and aligning, by a second subset of the hardwireddigital logic circuits of the integrated circuit, the plurality of readsto the one or more segments of the one or more genetic referencesequences.

In various embodiments, the method may additionally include accessing,by the integrated circuit on one or more of the plurality of physicalelectrical interconnects from a memory, the aligned plurality of reads.In such an instance the method may include sorting, by a third subset ofthe hardwired digital logic circuits of the integrated circuit, thealigned plurality of reads according to their positions in the one ormore genetic reference sequences. In certain instances, the method mayfurther include outputting, such as on one or more of the plurality ofphysical electrical interconnects of the integrated circuit, result datafrom the mapping and/or the aligning and/or the sorting, such as wherethe result data includes positions of the mapped and/or aligned and/orsorted plurality of reads.

Hence, in various instances, implementations of various aspects of thedisclosure may include, but are not limited to: apparatuses, systems,and methods including one or more features as described in detailherein, as well as articles that comprise a tangibly embodiedmachine-readable medium operable to cause one or more machines (e.g.,computers, etc.) to result in operations described herein. Similarly,computer systems are also described that may include one or moreprocessors and one or more memories coupled to the one or moreprocessors. Accordingly, computer implemented methods consistent withone or more implementations of the current subject matter can beimplemented by one or more data processors residing in a singlecomputing system or multiple computing systems. Such multiple computingsystems can be connected and can exchange data and/or commands or otherinstructions or the like via one or more connections, including but notlimited to a connection over a network (e.g. the Internet, a wirelesswide area network, a local area network, a wide area network, a wirednetwork, or the like), via a direct connection between one or more ofthe multiple computing systems, etc. A memory, which can include acomputer-readable storage medium, may include, encode, store, or thelike one or more programs that cause one or more processors to performone or more of the operations described herein.

The details of one or more variations of the subject matter describedherein are set forth in the accompanying drawings and the descriptionbelow. Other features and advantages of the subject matter describedherein will be apparent from the description and drawings, and from theclaims. While certain features of the currently disclosed subject matterare described for illustrative purposes in relation to an enterpriseresource software system or other business software solution orarchitecture, it should be readily understood that such features are notintended to be limiting. The claims that follow this disclosure areintended to define the scope of the protected subject matter.

DESCRIPTION OF DRAWINGS

The accompanying drawings, which are incorporated in and constitute apart of this specification, show certain aspects of the subject matterdisclosed herein and, together with the description, help explain someof the principles associated with the disclosed implementations. In thedrawings,

FIG. 1 is a block diagram of a hardware processor architecture inaccordance with an implementation.

FIG. 2 is a block diagram of a hardware processor architecture inaccordance with another implementation.

FIG. 3 is a block diagram of a hardware processor architecture inaccordance with yet another implementation

FIG. 4 shows a genetic sequence analysis pipeline.

FIG. 5 illustrates processing steps using a genetic sequence analysishardware platform.

FIG. 6 illustrates an apparatus in accordance with an implementation.

FIG. 7 illustrates an apparatus in accordance with an alternativeimplementation.

FIG. 8 illustrates a genomics processing system in accordance with animplementation.

When practical, similar reference numbers denote similar structures,features, or elements.

DETAILED DESCRIPTION

To address these and potentially other issues with currently availablesolutions, methods, systems, articles of manufacture, and the likeconsistent with one or more implementations of the current subjectmatter can, among other possible advantages, provide a sequence analysisapparatus for executing a sequence analysis pipeline on genetic sequencedata.

The following provides details of various implementations of a sequenceanalysis pipeline and platform.

In its most basic form, the body is comprised of cells, the cells formtissues, tissues form organs, organs form systems, and these systemsfunction together to ensure the body operates to sustain the life of theindividual. The cells of the body, therefore, are the building blocks oflife. More particularly, each cell has a nucleus, and within the nucleusof every cell reside chromosomes. Chromosomes are formed fromDeoxyribonucleic Acid, which has an organized but winding double helixstructure. The DNA itself is comprised of two opposed, but complementarystrands of nucleotides, which nucleotides comprise the genes that codefor the proteins that give the cells their structures and mediate thefunctions and regulations of the body's tissues and organs. Basically,proteins do most of the work of cells in maintaining the body's normalprocesses and functions.

Given the multiplicity of components of the body and the complexityinvolved in how they interact with one another to maintain the body'svarious processes and functions, there are a multiplicity of ways thatthe body may malfunction on any one of these different levels. Forinstance, in one such instance, there may be a malfunction in the way aparticular gene codes for a given protein, which dependent on theprotein and the nature of its malfunctioning can result in the onset ofa diseased state.

Accordingly, in diagnosing, preventing, and/or curing such diseasedstates, determining the genetic makeup of a subject may be extremelyuseful. For instance, once known, a person's genetic makeup, e.g., hisor her genomic composition, can be used for purposes of diagnosticsand/or for determining whether a person has or has the potential for adiseased state. Likewise, the knowledge of a person's genome may beuseful in determining various potential therapeutic modalities, such asdrugs, that can or cannot be used in a prophylactic or therapeuticregimen without causing harm to the user. In various instances,knowledge of a person's genome may also be employed to determine drugefficacy and/or problematic side effects of such drug use may bepredicted and/or identified. Potentially, the knowledge of a person'sgenome can be used to produce designer drugs, such as drugs tailor madeand optimized in accordance with a person's specific genetic makeup. Inparticular, in one instance, an engineered protein or nucleotidesequence can be fabricated to an individual's unique geneticcharacteristics so as to turn off or turn on the transcription of genesthat either over or under produce proteins and thereby amelioratediseased states.

Hence, in some instances, it is a goal of bioinformatics processing todetermine individual genomes of people, which determinations may be usedin gene discovery protocols as well as for prophylaxis and/ortherapeutic regimes to better enhance the livelihood of each particularperson and human kind as a whole. Further, knowledge of an individual'sgenome may be used such as in drug discovery and/or FDA trials to betterpredict with particularity which, if any, drugs will be likely to workon an individual and/or which would be likely to have deleterious sideeffects, such as by analyzing the individual's genome and/or a proteinprofile derived therefrom and comparing the same with predictedbiological response from such drug administration.

Such bioinformatics processing usually involves three well defined, buttypically separate phases of information processing. The first phaseinvolves DNA sequencing, where a subject's DNA is obtained and subjectedto various processes whereby the subject's genetic code is converted toa machine-readable digital code, e.g., a FASTQ file. The second phaseinvolves using the subject's generated digital genetic code for thedetermination of the individual's genetic makeup, e.g., determining theindividual's genomic nucleotide sequence. And the third phase involvesperforming one or more analyses on the subject's genetic makeup so as todetermine therapeutically useful information therefrom.

Preliminarily to Phase I, or primary processing, the genetic materialmust be preprocessed, so as to derive usable genetic sequence data. Thispreprocessing may be done manually or via an automated sequencer.Typically, preprocessing involves obtaining a biological sample from asubject, such as through venipuncture, hair, etc. and treating thesample to isolate the DNA therefrom. Once isolated the DNA may bedenatured, strand separated, and/or portions of the DNA may then bemultiplied, e.g., via polymerase chain reaction (PCR), so as to build alibrary of replicated strands that are now ready to be read, such as byan automated sequencer, which sequencer is configured to read thereplicate strands, e.g., by synthesis, and thereby determine thenucleotide sequences that makes up the DNA. Further, in variousinstances, such as in building the library of replicated strands, it maybe useful to provide for over-coverage when preprocessing a givenportion of the DNA. To perform this over-coverage, e.g., using PCR, mayrequire increased sample preparation resources and time, and thereforebe more expensive, but it often gives an enhanced probability of the endresult being more accurate.

Once the library of replicated strands has been generated they may beinjected into an automated sequencer that may then read the strands,such as by synthesis, so as to determine the nucleotide sequencesthereof. For instance, the replicated single stranded DNA may beattached to a glass bead and inserted into a test vessel, e.g., anarray. All the necessary components for replicating its complementarystrand, including labeled nucleotides, are also added to the vessel butin a sequential fashion. For example, all labeled “A”, “C”, “G”, and“T's” are added, either one at a time or all together to see which ofthe nucleotides is going to bind at position one. After each addition alight, e.g., a laser, is shone on the array. If the compositionfluoresces then an image is produced indicating which nucleotide boundto the subject location. More particularly, where the nucleotides areadded one at a time, if a binding event occurs, then its indicativefluorescence will be observed. If a binding event does not occur, thetest vessel may be washed and the procedure repeated until theappropriate one of the four nucleotides binds to its complement at thesubject location, and its indicative fluorescence is observed. Where allfour nucleotides are added at the same time, each may be labeled with adifferent fluorescent indicator, and the nucleotide that binds to itscomplement at the subject position may be determined, such as by thecolor of its fluorescence. This greatly accelerates the synthesisprocess.

Once a binding event has occurred, the complex is then washed and thesynthesis steps are repeated for position two. For example, a markednucleotide “A” may be added to the mix to determine if the complement atposition one is an “A”, and if so, all the sequences having thatcomplement will bind to the labeled “A” and will therefore fluoresce,and the samples will all be washed. Where the binding happened the boundnucleotide is not washed away, and then this will be repeated for allnucleotides for all positions until all the over-sampled nucleic acidsegments, e.g., reads, have been sequenced and the data collected.Alternatively, where all four nucleotides are added at the same time,each labeled with a different fluorescent indicator, only one nucleotidewill bind to its complement at the subject position, and the others willbe washed away, such that after the vessel has been washed, a laser maybe shone on the vessel and which nucleotide bound to its complement maybe determined, such as by the color of its fluorescence.

This continues until the entire strand has been replicated in thevessel. Usually a typical length of a sequence replicated in this manneris from about 100 to about 500 base pairs, such as between 150 to about400 base pairs, including from about 200 to about 350 base pairs, suchas about 250 base pairs to about 300 base pairs dependent on thesequencing protocol being employed. Further, the length of thesesegments may be predetermined, e.g., engineered, to accord with anyparticular sequencing machinery and/or protocol by which it is run. Theend result is a readout, or read, that is comprised of a replicated DNAsegment, e.g., from about 100 to about 1,000 nucleotides in length, thathas been labeled in such a manner that every nucleotide in the sequence,e.g., read, is known because of its label. Hence, since the human genomeis comprised of about 3.2 billion base pairs, and various knownsequencing protocols usually result in labeled replicated sequences,e.g., reads, from about 100 or 101 bases to about 250 or about 300 orabout 400 bases, the total amount of segments that need to be sequenced,and consequently the total number of reads generated, can be anywherefrom about 10,000,000 to about 40,000,000, such as about 15,000,000 toabout 30,000,000, dependent on how long the label replicated sequencesare. Therefore, the sequencer may typically generate about 30,000,000reads, such as where the read length is 100 nucleotides in length, so asto cover the genome once.

However, as indicated above, in such procedures, it may be useful tooversample the DNA such by about 5×, or about 10×, or about 20×, orabout 25×, or about 30×, or about 40×, or about 50×, or about 100×, orabout 200×, or about 250×, or about 500×, or about 1,000×, or about5,000×, or even about 10,000× or more, and as such the amount of primaryprocessing needed to be done and the time taken to do this can be quiteextensive. For instance, with 40X oversampling, wherein the varioussynthesized reads are designed to overlap to some extent, up to about1.2 billion reads may need to be synthesized. Typically, a largemajority if not all of these labeled sequences can be generated inparallel. The end result is that the initial biological genetic materialis processed, e.g., by sequencing protocols such as those summarizedherein, and a digital representation of that data is generated, whichdigital representation of data may be subjected to a primary processingprotocol. Particularly, the genetic material of a subject may bereplicated and sequenced in such a manner that a measurable electrical,radioactive and/or optical signal is generated, which signal is thenconverted, e.g., by the sequencer, into a digital representation of thesubject's genetic code. More particularly, primary processing mayinclude the conversion of images, such as recorded flashes of light orother electrical signal data, into FASTQ file data. Accordingly, thisinformation is stored as a FASTQ file, which may then be sent forfurther, e.g., secondary processing. A typical FASTQ file includes alarge collection of reads representing digitally encoded nucleotidesequences wherein each predicted base in the sequence has been calledand given a probability score that the called base at the indicatedposition is incorrect.

In many instances, it may be useful to further process the digitallyencoded sequence data obtained from the sequencer and/or sequencingprotocol, such as by subjecting the digitally represented data tosecondary processing. This secondary processing, for instance, can beused to assemble an entire genomic profile of an individual, such aswhere the individual's entire genetic makeup is determined, forinstance, where each and every nucleotide of each and every chromosomeis determined in sequential order such that the composition of theindividual's entire genome has been identified. In such processing, thegenome of the individual may be assembled such as by comparison to areference genome, such as a standard, e.g., one or more genomes obtainedfrom the human genome project, so as to determine how the individual'sgenetic makeup differs from that of the referent(s). This process iscommonly known as variant calling. As the difference between the DNA ofany one person to another is 1 in 1,000 base pairs, such a variantcalling process can be very labor and time intensive.

Accordingly, in a typical secondary processing protocol, a subject'sgenetic makeup is assembled by comparison to a reference genome. Thiscomparison involves the reconstruction of the individual's genome frommillions upon millions of short read sequences and/or the comparison ofthe whole of the individual's DNA to an exemplary DNA sequence model. Ina typical secondary processing protocol a FASTQ file is received fromthe sequencer containing the raw sequenced read data. For instance, incertain instances, there can be up to 30,000,000 reads or more coveringthe subject's genome, assuming no oversampling, such as where each readis about 100 nucleotides in length. Hence, in such an instance, in orderto compare the subject's genome to that of the standard referencegenome, it needs to be determined where each of these reads map to thereference genome, such as how each is aligned with respect to oneanother, and/or how each read can also be sorted by chromosome order soas to determine at what position and in which chromosome each readbelongs. One or more of these functions may take place prior toperforming a variant call function on the entire full-length sequence.Once it is determined where in the genome each read belongs, the fulllength genetic sequence may be determined, and then the differencesbetween the subject's genetic code and that of the referent can beassessed.

As the human genome is over 3 billion base pairs in length, efficientautomated sequencing protocols and machinery have been developed so asto effectuate the sequencing of such a genome within a time period thatcould be clinically useful. Such innovations in automated sequencinghave resulted in the capabilities of sequencing an entire genome in amatter of hours to days dependent on the number of genomes beingsequenced, the amount of oversampling involved, and the number ofprocessing resources being dedicated to the job. Hence, given theseadvancements in sequencing, a large amount of sequencing data is capableof being generated in a relatively short period of time. A result ofthese advancements, however, is the development of a bottleneck at thesecondary processing stage. In efforts to help overcome this bottleneckvarious software based algorithms have been developed to help expeditethe process of assembling a subject's sequenced DNA such as by areference based assembly process.

For instance, reference based assembly is a typical secondary processingassembly protocol involving the comparison of sequenced genomic DNA of asubject to that of one or more standards, e.g., known referencesequences. Various algorithms have been developed to help expedite thisprocess. These algorithms typically include some variation of one ormore of: mapping, aligning, and/or sorting the millions of readsreceived from the FASTQ file communicated by the sequencer, to determinewhere on each chromosome each particular read is located. Often a commonfeature behind the functioning of these various algorithms is their useof an index and/or an array to expedite their processing function.

For instance, with respect to mapping, a large quantity, e.g., all, ofthe sequenced reads may be processed to determine the possible locationsin the reference genome to which those reads could possibly align. Onemethodology that can be used for this purpose is to do a directcomparison of the read to the reference genome so as to find all thepositions of matching. Another methodology is to employ a prefix orsuffix array, or to build out a prefix or suffix tree, for the purposeof mapping the reads to various positions in the reference genome. Atypical algorithm useful in performing such a function is aBurrows-Wheeler transform, which is used to map a selection of reads toa reference using a compression formula that compresses repeatingsequences of data. A further methodology is to employ a hash table, suchas where a selected subset of the reads, a k-mer of a selected length“k”, e.g., a seed, are placed in a hash table as keys and the referencesequence is broken into equivalent k-mer portions and those portions andtheir location are inserted by an algorithm into the hash table at thoselocations in the table to which they map according to a hashingfunction. A typical algorithm for performing this function is “BLAST”, aBasic Local Alignment Search Tool. Such hash table based programscompare query nucleotide or protein sequences to one or more standardreference sequence databases and calculates the statistical significanceof matches. In such manners as these, it may be determined where anygiven read is possibly located with respect to a reference genome. Thesealgorithms are useful because they require less memory, fewer look ups,and therefore require fewer processing resources and time in theperformance of their functions, than would otherwise be the case, suchas if the subject's genome were being assembled by direct comparison,such as without the use of these algorithms.

Additionally, an aligning function may be performed to determine out ofall the possible locations a given read may map to on a genome, such asin those instances where a read may map to multiple positions in thegenome, which is in fact the location to which it actually was derived,such as by being sequenced therefrom by the original sequencingprotocol. This function may be performed on a number of the reads of thegenome and a string of ordered nucleotide bases representing a portionor the entire genetic sequence of the subject's DNA may be obtained.Along with the ordered genetic sequence a score may be given for eachnucleotide position, representing the likelihood that for any givennucleotide position, the nucleotide, e.g., “A”, “C”, “G”, “T” (or “U”),predicted to be in that position is in fact the nucleotide that belongsin that assigned position. Typical algorithms for performing alignmentfunctions are Needleman-Wunsch and Smith-Waterman. In either case, thesealgorithms perform sequence alignments between a string of the subject'squery genomic sequence and a string of the reference genomic sequencewhereby instead of comparing the entire genomic sequences, one with theother, segments of a selection of possible lengths are compared.

Once the reads have been assigned a position, such as relative to thereference genome, which may include identifying to which chromosome theread belongs and/or its offset from the beginning of that chromosome,the reads may be sorted by position. This may enable downstream analysesto take advantage of the oversampling described above. All of the readsthat overlap a given position in the genome will be adjacent to eachother after sorting and they can be organized into a pileup and readilyexamined to determine if the majority of them agree with the referencevalue or not. If they do not, a variant can be flagged.

Although these algorithms and the others like them go a ways toresolving the bottlenecks inherent in secondary processing, fasterperformance time and better accuracy are still desirable. Moreparticularly, although there has been advancement in the generation ofraw data, such as sequence data, the advancements in informationtechnologies have not kept up pace, leading to a data analysisbottleneck. This bottleneck is somewhat lessened by the development ofvarious algorithms, such as those described above, which help acceleratethese analyses, but there still exists a need for new technologies tohandle the computation, storage, and/or analysis of such data,especially as it relates to genomic sequence analysis, such as in asecondary processing stage.

For instance, employing standard protocols for performing secondaryprocessing on obtained genomic sequencing data, can take up to three (3)days or even up to a week or more to process the sequenced data so as togenerate clinically relevant genomic sequence information of anindividual. Employing various different optimized algorithms, such asthose described above, the time expended for secondary processing can bebrought down to a mere 27 to 48 hours. However, in order to achieve suchrapid results typically requires virtually all the generated reads,e.g., 30 million reads of 100 nucleotides each, to be processed inparallel and at the same time. Such parallel processing requiresextensive processing power involving massive CPU resources and stilltakes a relatively long time.

Further, in various instances, enhanced accuracy of results is desired.Such enhanced accuracy can be achieved through providing some amount ofoversampling of the sequenced genome. For example, as described above,it may be desirable to process the subject's DNA in such a manner thatat any given location of a sequence of nucleotides, there is anoversampling of that region. As indicated above, it may be desired tooversample any given region of the genome up to 10×, or 15×, or 20×, or25×, or 30×, or 40×, 50×, 100×, 250× or even 500× or 1,000 times ormore. However, where the genome is oversampled, such as by 40×, theamount of reads to be processed is roughly 30 Million×40 (dependent onthe length of the reads), which amounts to about 1.2 billion reads thatneed to be processed, when the entire genome is oversampled by 40×.Hence, although such oversampling typically results in greater accuracy,it is at a cost of taking more time and requiring more extensiveprocessing resources as each section of the genome is covered byanywhere from 1 to 40 times. Moreover, for certain oncology applicationsin which a clinician is trying to distinguish between the mutated genomeof cancer cells in the blood stream as distinct from the genome ofhealthy cells, oversampling of as much as 500×, or 1,000×, or 5,000×, oreven 10,000× may be employed.

The present disclosure, therefore, is directed to such new technologiesthat may be implemented in one or a series of genomics and/orbioinformatics protocols for performing genetic analysis, such assecondary processing, on obtained genomic sequencing data or a portionthereof. The sequencing data may be obtained directly from an automatedhigh throughput sequencer system, such as by a “Sequencing by Synthesis”454 automated sequencer from ROCHE, a HiSeq×Ten or a Solexia automatedsequencer from ILLUMINA, a “Sequencing by Oligonucleotide Ligation andDetection” (SOLiD) or Ion Torrent sequencer by LIFE TECHNOLOGIES, and/ora “Single Molecule Fluorescent Sequencing” sequencer by HELICOS GENETICANALYSIS SYSTEMS, or the like, such as by a direct linkage with thesequencing processing unit, or the sequencing data may be obtainedremotely, such as from a database, for instance, accessible via theinternet or other remote location accessible through a wirelesscommunications protocol, such as Wi-Fi, Bluetooth, or the like.

In certain aspects, these genetic analysis technologies may employimproved algorithms that may be implemented by software that is run in aless processing intensive and/or less time consuming manner and/or withgreater percentage accuracy. For instance, in certain embodiments,improved algorithms for performing such secondary processing, asdisclosed herein, is provided. In various particular embodiments, theimproved algorithms are directed to more efficiently and/or moreaccurately performing one or more of mapping, aligning, and/or sortingfunctions, such as on a digital representation of DNA sequence dataobtained from a sequencing platform, such as in a FASTQ file formatobtained from an automated sequencer such as one of those set forthabove.

In certain embodiments, improved algorithms directed to more efficientlyand/or more accurately performing one or more of local realignment,duplicate marking, base quality score recalibration, variant calling,compression, and/or decompression functions are provided. Further, asdescribed in greater detail herein below, in certain aspects, thesegenetic analysis technologies may employ on or more algorithms, such asimproved algorithms, that may be implemented by hardware that is run ina less processing intensive and/or less time consuming manner and/orwith greater percentage accuracy than various software implementationsfor doing the same.

In particular embodiments, a platform of technologies for performinggenetic analyses are provided where the platform may include theperformance of one or more of: mapping, aligning, sorting, localrealignment, duplicate marking, base quality score recalibration,variant calling, compression, and/or decompression functions. In certaininstances, the implementation of one or more of these platform functionsis for the purpose of performing one or more of determining and/orreconstructing a subject's consensus genomic sequence, comparing asubject's genomic sequence to a referent sequence, e.g., a reference ormodel genetic sequence, determining the manner in which the subject'sgenomic DNA differs from a referent, e.g., variant calling, and/or forperforming a tertiary analysis on the subject's genomic sequence, suchas for genome-wide variation analysis, gene function analysis, proteinfunction analysis, e.g., protein binding analysis, quantitative and/orassembly analysis of genomes and/or transcriptomes, as well as forvarious diagnostic, and/or a prophylactic and/or therapeutic evaluationanalyses.

Further, in various embodiments, a bioinformatics processing regime, asdisclosed herein, may be employed for the purpose of creating one ormore masks, such as a genome reference mask, a default mask, a diseasemask, and/or an iterative feed back mask, which may be added to themapper and/or aligner, e.g., along with a reference, wherein the maskset is configured so as to identify a particular area or object ofinterest. For instance, in one embodiment, the methods and apparatusesdescribed herein may be employed so as to create genome reference mask,such as by creating a mask-set that can be loaded into the mapper and/oraligner along with a reference, wherein the mask set is configured so asto identify areas of high importance and/or relevance, e.g., to thepractitioner or subject, and/or so as to identify areas having increasedsusceptibility to errors. In various embodiments, the mask-set mayprovide intelligent guidance to the mapper and/or aligner such as onwhich areas of the genome to focus on to improve quality. Masks,therefore, can be created in a layered manner to provide varying levelsor iterations of guidance based on various specific applications. Eachmask accordingly could identify the areas of interest and provide aminimum quality target for the area. Additionally, a default mask may beemployed to provide guidance, such as on an identified, e.g., typical,“high value” areas of the genome. Such areas could include known codingareas, control areas, etc. as well as areas that are well known toproduce errors. Further, a disease mask, or application specific mask,may be employed to the mask-set that identifies areas of highimportance, such as areas that require very high levels of accuracybased on known markers, e.g., Cancer. Further still, iterative feedbackmasking may be employed, such as by adding a new, ad-hoc mask, that maybe specifically designed by using feedback from a tertiary analysissystem (like Cypher Genomics) that has identified areas of concern basedon observed errors or inconsistencies.

As indicated above, in one aspect one or more of these platformfunctions, e.g., mapping, aligning, sorting, realignment, duplicatemarking, base quality score recalibration, variant calling, compression,and/or decompression functions is configured for implementation insoftware. In another embodiment, one or more of these platformfunctions, e.g., mapping, aligning, sorting, local realignment,duplicate marking, base quality score recalibration, decompression,variant calling, compression, and/or decompression functions isconfigured for implementation in hardware.

Accordingly, in certain instances, methods are presented herein wherethe method involves the performance of an algorithm, such as analgorithm for performing one or more genetic analysis functions such asmapping, aligning, sorting, realignment, duplicate marking, base qualityscore recalibration, variant calling, compression, and/or decompressionwhere the algorithm has been optimized in accordance with the manner inwhich it is to be implemented. In particular, where the algorithm is tobe implemented in a software solution, the algorithm and/or itsattendant processes, has been optimized so as to be performed fasterand/or with better accuracy for execution by that media. Likewise, wherethe functions of the algorithm are to be implemented in a hardwaresolution, the hardware has been designed to perform these functionsand/or their attendant processes in an optimized manner so as to beperformed faster and/or with better accuracy for execution by thatmedia. These methods, for instance, can be employed such as in aniterative variant calling procedure.

Hence, in one aspect, presented herein are systems, apparatuses, andmethods for implementing bioinformatic protocols, such as for performingone or more functions for analyzing genetic data, such as genomic data,for instance, via one or more optimized algorithms and/or on one or moreoptimized integrated circuits, such as on one or more hardwareprocessing platforms. Hence, in one instance, systems and methods areprovided for implementing one or more algorithms for the performance ofone or more steps for analyzing genomic data in a bioinformaticsprotocol, such as where the steps may include the performance of one ormore of: mapping, aligning, sorting, local realignment, duplicatemarking, base quality score recalibration, variant calling, compression,and/or decompression. In another instance, systems and methods areprovided for implementing the functions of one or more algorithms forthe performance of one or more steps for analyzing genomic data in abioinformatics protocol, as set forth herein, wherein the functions areimplemented on a hardware accelerator, which may or may not be coupledwith one or more general purpose processors and/or super computers.

More specifically, in some instances, methods for performing secondaryanalytics on data pertaining to the genetic composition of a subject areprovided. In one instance, the analytics to be performed may involvereference based reconstruction of the subject genome. For instance,referenced based mapping involves the use of a reference genome, whichmay be generated from sequencing the genome of a single or multipleindividuals, or it may be an amalgamation of various people's DNA thathave been combined in such a manner so as to produce a prototypical,standard reference genome to which any individual's DNA may be compared,for example, so as to determine and reconstruct the individual's geneticsequence and/or for determining the difference between their geneticmakeup and that of the standard reference, e.g., variant calling.

More particularly, a reason for performing a secondary analysis on asubject's sequenced DNA is to determine how the subject's DNA variesfrom that of the reference. More specifically, to determine one, amultiplicity, or all the differences in the nucleotide sequence of thesubject from that of the reference. For instance, the differencesbetween the genetic sequences of any two random persons is 1 in 1,000base pairs, which when taken in view of the entire genome of over 3billion base pairs amounts to a variation of up to 3,000,000 divergentbase pairs per person. Determining these differences may be useful suchas in a tertiary analysis protocol, for instance, so as to predict thepotential for the occurrence of a diseased state, such as because of agenetic abnormality, and/or the likelihood of success of a prophylacticor therapeutic modality, such as based on how a prophylactic ortherapeutic is expected to interact with the subject's DNA or theproteins generated therefrom. In various instances, it may be useful toperform both a de novo and a reference based reconstruction of thesubject's genome so as to confirm the results of one against the other,and to, where desirable, enhance the accuracy of a variant callingprotocol.

In various instances, as set forth above, it may be useful in performinga primary sequencing protocol to produce oversampling for one or moreregions of the subject's genome. These regions may be selected based onknown areas of increased variability, suspected regions of variability,such as based on the condition of the subject, and/or on the entiregenome generally. In its basic form, as indicated above, based on thetype of sequencing protocols performed, sequencing produces readouts,e.g., reads, that are digital representations of the subject's geneticsequence code. These read lengths are typically designed based on thetype of sequencing machinery being employed. For instance, the 454automated sequencer from ROCHE, typically produces read lengths from 100or 150 base pairs in length to about 1,000 base pairs; for ILLUMINA theread lengths are typically engineered to be from about 100 or 101 toabout 150 base pairs in length for some of their technology, and 250base pairs in length for other of their technology; for LIFETECHNOLOGIES the read lengths are typically engineered to be from about50 to about 60 base pairs in length for their SOLiD technology and from35 to 450 base pairs in length for their Ion Torrent technology; and forthe HELICOS GENETIC ANALYSIS SYSTEMS the read lengths may vary but maytypically be less than 1,000 nucleotides in length.

However, because the processing of the DNA sample required to produceengineered read lengths of a specific size is both labor and chemistryintensive, and because the sequencing itself often depends on thefunctioning of the sequencing machinery, there is some possibility thaterrors may be made throughout the sequencing process thereby introducingan abnormality into that portion of the sequenced genome where the erroroccurred. Such errors can be problematic especially where a purpose forreconstructing the subject's genome is to determine how it or at least aportion of the genome varies from a standard or model reference. Forinstance, a machine or chemistry error resulting in the change of onenucleotide, e.g., in a read, for another will give a false indication ofa variation that is not really there. This can result in an incorrectvariant call and may further result in the false indication of adiseased state and the like. Accordingly, because of the possibility ofmachine, chemistry, and/or even human error in the execution of asequencing protocol, in many instances, it is desirable to buildredundancy into an analysis system, such as by oversampling portions ofor the entire genome. More particularly, as an automated sequencerproduces a FASTQ file calling out a sequence of reads having nucleotidesat a given position along with the probability that the call for a givennucleotide being at the called position is actually incorrect, e.g., abase call, it is often desirable to employ methods, such asoversampling, for ensuring that base calls made by the sequencingprocesses can be detected and corrected.

Hence, in performing the methods herein described, in certain instances,a primary sequencing protocol is performed in such a manner so as toproduce a sequenced genome where a portion or the entire genome isoversampled by about 10×, about 15×, about 20×, about 25×, about 30×about 40×, such as about 50× or more. Accordingly, where the readlengths are engineered to be about 50-60 base pairs in length, thisoversampling can result in about 2 to about 2.5 billion reads, or wherethe read lengths are about 100 or 101 base pairs in length, oversamplingmay result in about 1 to about 1.2 billion reads, and where the readlengths are about 1,000 base pairs in length, about 50 to about 100million reads may be generated by the sequencer, such as where theoversampling is about 40×. More particularly, in such an instance,because of the 40× oversampling, at any given point in the genome it isexpected that there will be 40 reads to cover any one position albeit,the given position might be at the beginning of one read, the middle ofanother, and the end of another, but it is expected to be covered about40 times.

Therefore, such oversampling produces regions of the sequenced genomethat are covered by a multiplicity of reads, e.g., duplications, such asup to about 40 reads, for instance, where the oversampling is about 40×.These at least partial duplications are useful in determining whetherany given variation in any particular read is in fact an actual genomicvariation or rather a machine or chemistry artifact. Hence, oversamplingcan be employed to improve the accuracy in reconstructing the subject'sgenome, especially in instances where the subject's genome is to becompared against a reference genome so as to determine those instanceswhere the subject's genetic sequence differs from that of the referencegenetic sequence. In a manner such as this, as described in greaterdetail herein below, it can be confirmed that any given variationbetween the reconstructed sequence and the model is in fact due to thepresence of an actual variant and not an error in the initial processingof sample DNA, or read alignment software, etc.

For instance, in building the genetic sequence of the individual'ssequenced DNA, it must be determined what nucleotide goes where in thegrowing string of nucleotides. In order to determine what nucleotidegoes where, the various reads can be organized and a pile up of readscovering duplicate locations can be built up. This allows for acomparison to be made of all the reads covering the same locations so asto more accurately determine if there is an actual variation at anygiven position or if there may be an error in any one read at theposition in question in the pileup. For example, if there is only one ortwo of the reads out of the 40 that has a particular nucleotide atposition X, and all 38 or 39 other reads agree on a different nucleotidebeing at that position, then the two outlying reads may be excluded asbeing in error, at least at this specific location.

More particularly, where there are a multiplicity of reads generated forany one location of the subject's genome, there are likely to bemultiple overlaps or pile-ups for any given nucleotide position. Thesepile-ups represent the coverage for any particular location and may beuseful for determining with better accuracy the correct sequence of thesubject's genome. For instance, as indicated, sequencing results in theproduction of reads, and in various instances, the reads produced areover sampled, and so at various positions various particular reads willoverlap. This overlapping is useful for determining the actual samplegenome such as with a high probability of correctness.

The purpose, therefore, may be to scan over the reference genomeincrementally multiple times, as described in greater detail hereinbelow, so as to more accurately reconstruct the subject's genome, andwhere it is desirable to determine how the subject's genome differs froma different genome, e.g., a model genome, the use of pile-ups can moreaccurately identify errors, such as chemical, machine, or read errors,and distinguish them from actual variants. More specifically, where thesubject has an actual variation at position X, the majority of reads inthe pile up should verify, e.g., include, that variation. Statisticalanalysis procedures, such as those described herein, may then performedto determine the actual genetic sequence of the subject with all itsvariants from a reference genome.

For instance, where the subject's genetic sequence is to be rebuilt withrespect to the use of a reference genome, once the reads, e.g., apile-up of reads, have been generated, the next steps may be to mapand/or align and/or sort the reads to one or more reference genomes(e.g., the more exemplary reference genomes available as models thebetter the analysis is likely to be) and thereby rebuild the genome ofthe subject, this results in a series of reads that have been mappedand/or aligned with the reference genome(s) at all possible positionsalong the chain where there is a match, and at each such position theyare given a probability score as to the probability that they actuallybelong in that position.

Accordingly, in various instances, once the reads have been generated,their positions mapped, e.g., the potential locations in the referencegenome to which the reads may map have been determined, and theirsequential order aligned, the actual genetic sequence of the subject'sgenome may be determined, such as by performing a sorting function onthe aligned data. Further, once the actual sample genome is known andcompared to the reference genome, the variations between the two can bedetermined, a list of all the variations/deviations between thereference genome and the sample genome are determined and called out.Such variations between the two genetic sequences may be due to a numberof reasons.

For instance, there may be a single nucleotide polymorphism (SNP), suchas wherein one base in the subject's genetic sequence has beensubstituted for another; there may be more extensive substitutions of aplurality of nucleotides; there may be an insertion or a deletion, suchas where one or a multiplicity of bases have been added to or deletedfrom the subject's genetic sequence, and/or there may be a structuralvariant, e.g., such as caused by the crossing of legs of twochromosomes, and/or there may simply be an offset causing a shift in thesequence. In various instances, a variant call file containing all thevariations of the subject's genetic sequence to the reference sequenceis generated. More particularly, in various embodiments, the methods ofthe disclosure include generating a variant call file (VCF) identifyingone or more, e.g., all of the genetic variants in the individual whoseDNA was sequenced, e.g., relevant to one or more reference genomes. TheVCF in its basic form is a list of locations of variants and their type:e.g., chromosome 3, at position X, an “A” is substituted for a “T”, etc.

However, as indicated above, in order to generate such a file, thegenome of the subject must be sequenced and rebuilt prior to determiningits variants. There are, however, several problems that may occur whenattempting to generate such an assembly. As noted above, there may beproblems with the chemistry, the sequencing machine, and/or human errorthat occurs in the sequencing process. Additionally, there may begenetic artifacts that make such reconstructions problematic. Forinstance, a problem with performing such assemblies is that there aresometimes huge portions of the genome that repeat themselves, such aslong sections of the genome that include the same strings ofnucleotides. Hence, because any genetic sequence is not uniqueeverywhere, it may be difficult to determine where in the genome anidentified read actually maps and aligns.

For instance, dependent on the sequencing protocol employed shorter orlonger reads may be produced. Longer reads are useful in that the longerthe read the less likely it is to show up in multiple locations in thegenome. Having fewer possible locations to evaluate can also speed upthe system. However, the longer the reads the more problematic they maybe because the more likely they are to include a real or falsevariation, e.g., caused by an SNP, InDel (insertion or deletion), or amachine error, or the like, resulting in a no match between the read andthe reference genome. On the other hand, shorter reads are usefulbecause the shorter the read the less likely it is to cover a positionthat codes for a variant. A problem with shorter reads however is thatthe shorter the read the more likely it is to show up at multiplepositions in the genome, thus requiring additional processing time andresources so as to determine which out of all possible positions is themost likely actual position to where it aligns. Ideally what may beachieved, such as by practicing the methods herein disclosed, is that avariant call file may be produced wherein a list of the sequenced genome(the query sequence) is generated that shows where all the variant basepairs are, making sure each variant called is an actual variant and notsimply a chemistry or machine read or other human based error.

There are, therefore, two main possibilities for variation. For one,there is an actual variation at the particular location in question, forinstance, where the person's genome is in fact different at a particularlocation than that of the reference, e.g., there is a natural variationdue to an SNP (one base substitution), an Insertion or Deletion (of oneor more nucleotides in length), and/or there is a structural variant,such as where the DNA material from one chromosome gets crossed onto adifferent chromosome or leg, or where a certain region gets copied twicein the DNA. Alternatively, a variation may be caused by there being aproblem in the read data, either through chemistry or the machine,sequencer or aligner, or other human error. Accordingly, the methodsdisclosed herein may be employed in a manner so as to compensate forthese types of errors, and more particularly so as to distinguish errorsin variation due to chemistry, machine or human, and real variations inthe sequenced genome. More specifically, the methods, apparatuses, andsystems for employing the same, as here in described, have beendeveloped so as to clearly distinguish between these two different typesof variations and therefore to better ensure the accuracy of any callfiles generated so as to correctly identify true variants.

Further, in various embodiments, once the subject's genome has beenreconstructed and/or a VCF has been generated, such data may then besubjected to tertiary processing so as to interpret it, such as fordetermining what the data means with respect to identifying whatdiseases this person may or may have the potential for suffer fromand/or for determining what treatments or lifestyle changes this subjectmay want to employ so as to ameliorate and/or prevent a diseased state.For example, the subject's genetic sequence and/or their variant callfile may be analyzed to determine clinically relevant genetic markersthat indicate the existence or potential for a diseased state and/or theefficacy of a proposed therapeutic or prophylactic regimen may have onthe subject. This data may then be used to provide the subject with oneor more therapeutic or prophylactic regimens so as to better thesubject's quality of life, such as treating and/or preventing a diseasedstate.

More particularly, medical science technologies have advanced inconjunction with the advancement of information technologies, whichadvancement has enhanced our ability to store and analyze medical data.Hence, once one or more of an individual's genetic variations aredetermined, such variant call file information can be used to developmedically useful information, which in turn can be used to determine,e.g., using various known statistical analysis models, health relateddata and/or medical useful information, e.g., for diagnostic purposes,e.g., diagnosing a disease or potential therefore, clinicalinterpretation (e.g., looking for markers that represent a diseasevariant), whether the subject should be included or excluded in variousclinical trials, and other such purposes. As there are a finite numberof diseased states that are caused by genetic malformations, in tertiaryprocessing variants of a certain type, e.g., those known to be relatedto the onset of diseased states, can be queried for, such as bydetermining if one or more genetic based diseased markers are includedin the variant call file of the subject.

Consequently, in various instances, the methods herein disclosed mayinvolve analyzing, e.g., scanning, the VCF and/or the generatedsequence, against a known disease sequence variant, such as in a database of genomic markers therefore, so as to identify the presence of thegenetic marker in the VCF and/or the generated sequence, and if presentto make a call as to the presence or potential for a genetically induceddiseased state. As there are a large number of known genetic variationsand a large number of individual's suffering from diseases caused bysuch variations, in some embodiments, the methods disclosed herein mayentail the generation of one or more databases linking sequenced datafor an entire genome and/or a variant call file pertaining thereto,e.g., such as from an individual or a plurality of individuals, and adiseased state and/or searching the generated databases to determine ifa particular subject has a genetic composition that would predisposethem to having such diseased state. Such searching may involve acomparison of one entire genome with one or more others, or a fragmentof a genome, such as a fragment containing only the variations, to oneor more fragments of one or more other genomes such as in a database ofreference genomes or fragments thereof.

Further, it is understood that the genetic sequences to be employed inthese manners may be DNA, ssDNA, RNA, mRNA, rRNA, tRNA, or the like.Hence, although throughout the present disclosure various mention ismade to various methods and apparatuses for analyzing genomic DNA, invarious instances, the systems, apparatuses and methods disclosed hereinare equally suitable for performing their respective functions, e.g.,analysis, on all types of genetic material including DNA, ssDNA, RNA,mRNA, rRNA, tRNA, and the like. Additionally, in various instances, themethods of the disclosure may include analyzing the generated geneticsequence, e.g., DNA, ssDNA, RNA, mRNA, rRNA, tRNA, and the like, fromthe subject and determining therefrom the protein variations which arelikely to be caused by the genetic sequence and/or determining and/orpredicting the potential for a diseased state therefrom, such as due toan error in protein expression. It is to be noted that the geneticsequence obtained can represent an intron or an exon, for instance, thegenetic sequence can be for a coding portion of the DNA only, such aswhere an exome is obtained and using known processing techniques onlythe coding regions, or non-coding regions, may be sequenced, which canlead to faster sequencing and/or faster processing times, albeitinvolving a more difficult sample preparation procedure.

Currently, such steps and analyses herein described are typicallyperformed in various distinct and unrelated steps often employingdifferent analytic machines at different locations. Accordingly, invarious aspects the methods and systems of the disclosure are performedby a single apparatus and/or at one location, such as in conjunctionwith an automated sequencer or other apparatus configured to generategenetic sequence data. In various instances, a plurality of apparatusesmay be employed at the same location, or a multiplicity of remotelocations, and in some instances, the methods may involve two or moreprocessing units being deployed at two or more locations.

For instance, in various aspects a pipeline may be provided wherein thepipeline includes performing one or more analytic functions, asdescribed herein, on a genomic genetic sequence of one or moreindividuals, such as data obtained in a digital, e.g., FASTQ, fileformat from an automated sequencer. A typical pipeline to be executedmay include one or more of sequencing genetic material, such as aportion or an entire genome, of one or more subjects, which geneticmaterial may include DNA, ssDNA, RNA, rRNA, tRNA, and the like, and/orin some instances the genetic material may represent coding ornon-coding regions, such as exomes, episomes of the DNA. The pipelinemay include one or more of performing a base calling and/or errorcorrection operation, such as on the digitized genetic data, and/or mayinclude one or more of performing a mapping, an alignment, and/or asorting function on the genetic data. In certain instances, the pipelinemay include performing one or more of a realignment, a deduplication, abase quality or score recalibration, a reduction and/or compression,and/or a decompression on the digitized genetic data. In certaininstances the pipeline may include performing a variant callingoperation on the genetic data.

Therefore, in various instances, a pipeline of the disclosure mayinclude one or more modules, wherein the modules are configured forperforming one or more functions, such as a base calling and/or errorcorrection operation and/or a mapping and/or an alignment and/or asorting function on genetic data, e.g., sequenced genetic data. And invarious instances, the pipeline may include one or more modules, whereinthe modules are configured for performing one more of a localrealignment, a deduplication, a base quality score recalibration, avariant calling, a reduction and/or compression, and/or a decompressionon the genetic data. Many of these modules may either be performed bysoftware or on hardware or remotely, e.g., via software or hardware,such as on the cloud or a remote server and/or server bank.

Additionally, many of these steps and/or modules of the pipeline areoptional and/or can be arranged in any logical order and/or omittedentirely. For instance, the software and/or hardware disclosed hereinmay or may not include a base calling or sequence correction algorithm,such as where there may be concern that such functions may result in astatistical bias. Consequently the system will either include or willnot include the base calling and/or sequence correction function,respectively, dependent on the level of accuracy and/or efficiencydesired. And as indicated above, one or more of the pipeline functionsmay be employed in the generation of a genomic sequence of a subjectsuch as through a reference based genomic reconstruction. Also asindicated above, in certain instances, the output from the pipeline is avariant call file indicating a portion or all the variants in a genomeor a portion thereof.

Accordingly, as indicated above, the output of performing a sequencingprotocol, such as one or more of those set forth above, is typically adigital representation of the subject's genetic material, such as in aFASTQ file format. However, an autorad that has been digitallytranscribed may also be employed. More particularly, the output from asequencing protocol may include a plurality of reads, where each readincludes a sequence, e.g., a string, of nucleotides where the positionof every nucleotide has been called, and a quality score representingthe probability that the called nucleotide is wrong. However, thequality of these outputs may be improved by various pre-processingprotocols so as to achieve higher quality of scores, which one or moreof such protocols may be employed in the methods disclosed herein.

For instance, in certain instances, the raw FASTQ file data may beprocessed to clean up the initial base calls obtained from thesequencer/reader, such as in a primary processing stage, e.g., prior tothe secondary processing described herein above. Specifically, thesequencer/reader typically analyzes the sequencing data, such as thefluorescent data indicating which nucleotide is at what position, andconverts the image data into a base call with a quality score, such aswhere the quality score is based on the comparative brightness of thefluorescence at each position. A specialized algorithm may be employed,such as in a primary processing stage, to correctly analyze thesedistinctions in fluorescence so as to more accurately make theappropriate base call. As indicated above, this step may be included ina pipeline of steps and may be implemented via software or hardware orboth, however, in this instance would be part of a primary processingplatform.

An additional preprocessing step may include an error correctionfunction, which may include an attempt to take the millions to billionsof reads in the FASTQ file and correct some proportion of any mechanicalsequencing error with the information pertaining to the base call andquality score available prior to any further processing such as mapping,alignment, and/or sorting functions, etc. For instance, the reads withinthe FASTQ file may be analyzed to determine if there are anysub-sequences in any of the reads that appear in other reads, whichbecause of the duplicate coverage can increase confidence that thesubsequences in the reads may be correct. This may be implemented bybuilding a hash table containing all possible k-mers of a selectedlength, k, from every read, and storing with each one its frequency andalso which bases immediately follow it and with what probability. Then,using the hash table each read can be rescanned. As each k-mer in aparticular read is looked up in the hash table, and evaluation can bemade as to whether the base immediately following that k-mer is likelyto be correct or not. If it is unlikely, then it can be replaced withthe most likely one to follow from the table. Subsequent k-mers for thatread will then include the corrected base as the value at that positionand the process is repeated. This can be highly effective in correctingerrors because oversampling enables gathering accurate statistics forpredicting what comes next after each k-mer. However, as indicatedabove, such corrections could add statistical biasing to the system,such as due to false corrections, to the data, and so these procedurescan be skipped if desired.

Accordingly, in accordance with the aspects of the disclosure, invarious instances, the methods, apparatuses, and/or systems of thedisclosure, may include obtaining read data, that either have or havenot been preprocessed, such as by being obtained directly from a FASTQfile of an automated sequencer, and subjecting the obtained data to oneor more of a mapping, aligning, and/or sorting function. The performanceof such functions may be useful, for instance, because, as set forthabove, in various instances, the sequencing data typically generated byvarious automated sequencers, e.g., reads, have lengths that aresubstantially shorter than the entire genomic sequence being analyzed,and since the human genome typically has a multiplicity of repetitivesections, and is known to have various repeating patterns in it, theremay be therefore a multiplicity of locations that any given readsequence may correspond to a segment in the human genome. Consequently,given all the possibilities a given read may match to the sequence ofthe genome, such as because of various repeating sequences in thegenome, etc. the raw read data may not clearly indicate which one of thepossibilities is in fact the correct location from which it was derived.Hence, for each read it will need to be determined to where in thegenome the reads actually map. Additionally, it may also be useful todetermine the sequential alignment of the reads, so as to determine theactual sequence identity of the subject, and/or it may also be useful todetermine the chromosomal location for each portion of the sequence.

Accordingly, in various instances, the methods of the disclosure may bedirected to mapping, aligning, and/or sorting the raw read data of theFASTQ files so as to find all the likely places that a given read may bealigned, and/or to determine the actual sequence identify of a subject,and/or to determine the chromosome location for each portion of thesequence. For example, mapping may be employed so as to map thegenerated reads to the reference genome and thereby find the locationwhere each read appears to match well to the genome, e.g., finding allthe places where there might be a good score for aligning any given readto the reference genome. Mapping therefore may involve taking one ormore, e.g., all, of the raw or preprocessed reads received from theFASTQ file and comparing the reads with one or more reference genomesand determining where the read may match with the reference genome(s).In its basic from, mapping involves finding the location(s) in thereference genome where one or more of the FASTQ reads obtained from thesequencer appears to match.

Likewise, alignment may be employed so as to evaluate all the candidatelocations of the individual reads against a window of the referencegenome to determine where and how the read sequences best align to thegenome. However, performing an alignment may be difficult due tosubstitutions, insertions, deletions, structural variations, and thelike which may prevent the read from aligning exactly. There are,therefore, several different ways to get an alignment, but to do so mayrequire making changes in the read, where each change that needs to bemade to get the appropriate alignment results in a lower confidencescore. For instance, any given read may have substitutions, insertions,and/or deletions as compared to the reference genome, and thesevariations need to be accounted for in performing an alignment.

Accordingly, along with the predicted alignment a probability score thatthe predicted alignment is correct may also be given. This scoreindicates the best alignment for any given read amongst multiplelocations where that read may align. For example, the alignment score ispredicated upon how well a given read matches a potential map locationand may include stretching, condensing, and changing bits and pieces ofthe read so as to get the best alignment.

The score will reflect all the ways the read was changed so as toaccommodate the reference. For instance, in order to generate analignment between the read and the reference one or more gaps in theread may need to be inserted, wherein the insertion of each gaprepresents a deletion in the read over the reference. Likewise,deletions may need to be made in the read, wherein each deletionrepresents an insertion in the read over the reference. Additionally,various bases may need to be changed such as due to one or moresubstitutions. Each of these changes are made to make the read(s) moreexactly align to the reference, but each change comes with a cost to thequality score, which score is a measure as to how well the entire readmatches to some region of the reference. The confidence in such qualityscores is then determined by looking at all the locations the read canbe made to map to the genome and comparing the scores at each location,and choosing the one with the highest score. More particularly, wherethere are multiple positions with high quality scores, then confidenceis low, but where the difference between the first and second bestscores is large, then confidence is high. At the end, all the proposedreads and confidence scores are evaluated and the best fit is selected.

Once the reads are assigned a position relative to the reference genome,which consists of identifying to which chromosome the read belongs andits offset from the beginning of that chromosome, they may be sorted,such as by position. This enables downstream analyses to take advantageof the various oversampling protocols described herein. All of the readsthat overlap a given position in the genome maybe be adjacent to eachother after sorting and they can be piled up and readily examined todetermine if the majority of them agree with the reference value or not.If they do not, as indicated above, a variant can be flagged.

As indicated above, the FASTQ file obtained from the sequencer iscomprised of a plurality, e.g., millions to a billion or more, of readsconsisting of short strings of nucleotide sequence data representing aportion or the entire genome of an individual. Mapping, in general,involves plotting the reads to all the locations in the reference genometo where there is a match. For example, dependent on the size of theread there may be one or a plurality of locations where the readsubstantially matches a corresponding sequence on the reference genome.Accordingly, the mapping and/or other functions disclosed herein may beconfigured for determining where out of all the possible locations oneor more reads may match to in the reference genome is actually the truelocation to where they map.

It is possible to compare every read with every position in the 3.2billion reference genome to determine where, if any, the reads match tothe reference genome. This may be done, for instance, where the readlengths approach about 100,000 nucleotides, about 200,000 nucleotides,about 400,000 nucleotides, about 500,000 nucleotides, even about1,000,000 or more nucleotides in length. However, where the reads aresubstantially shorter in length, such as where there are 50 millionreads or more, e.g., 1 billion reads, this process could take a verylong time and require a large amount of computing resources.Accordingly, there are several methods, such as described herein, thathave been developed for aligning the FASTQ reads to the reference genomein a much quicker manner. For instance, as disclosed above, one or morealgorithms may be employed so as to map one or more of the readsgenerated by the sequencer, e.g., in a FASTQ file, and match them to thereference genome, so as to determine where in the reference genome thesubject reads potentially map.

For instance, in various methods, an index of the reference isgenerated, so that the reads or portions of the reads may be looked upin the index, retrieving indications of locations in the reference, soas to map the reads to the reference. Such an index of the reference canbe constructed in various forms and queried in various manners. In somemethods, the index may include a prefix and/or a suffix tree. In othervarious methods, the index may include a Burrows/Wheeler transform ofthe reference. In further methods, the index may include one or morehash tables, and a hash function may be performed on one or moreportions of the reads in an effort to map the reads to the reference. Invarious instances, one or more of these algorithms may be performedsequentially or at the same time so as to accurately determine where oneor more, e.g., a substantial portion or every, read correctly matcheswith the reference genome.

Each of these algorithms may have advantages and/or disadvantages. Forexample, a prefix and/or suffix Tree and/or a Burrows/Wheelertransformation may be performed on the sequence data in such a mannerthat the index of the reference genome is constructed and/or queried asa tree-like data structure, where starting from a single-base or shortsubsequence of a read, the subsequence is incrementally extended withinthe read, each incremental extension stimulating accesses to the index,tracing a path through the tree-like data structure, until thesubsequence becomes unique enough, e.g., an optimal length has beenattained, and/or a leaf node is reached in the tree-like data structure,the leaf or last-accessed tree node indicating one or more positions inthe reference genome from which the read may have originated. Thesealgorithms, therefore, typically do not have a fixed length for the readsubsequences that may be mapped by querying the index. A hash function,however, often employs a fixed length comparison unit that may be theentire length of the read, but is often times a length that is somesub-portion thereof, which sub-portion is termed a seed. Such seeds canbe shorter or longer, but unlike with the prefix and/or suffix treesand/or the Burrows/Wheeler transformations, the seeds of the readsemployed in a hash function are typically of a preselected, fixedlength.

A prefix and/or suffix tree is a data structure that is built up fromthe reference genome, such that each link from a parent node to a childnode is labeled or associated with a nucleotide or sequence ofnucleotides, and each path from a root node through various links andnodes traces a path whose associated aggregate nucleotide sequencematches some continuous subsequence of the reference genome. The nodereached by such a path is implicitly associated with the referencesubsequence traced by its path from the root. Proceeding from the rootnode, subsequences in a prefix tree grow forward in the referencegenome, whereas subsequences in a suffix tree grow backward in thereference genome. Both a prefix tree and a suffix tree may be used in ahybrid prefix/suffix algorithm, so that subsequences may grow in eitherdirection. Prefix and suffix trees may also contain additional links,such as jumping from a node associated with one reference subsequence toanother node associated with a shorter reference subsequence.

For instance, a tree-like data structure serving as an index of thereference genome may be queried by tracing a path through the tree,corresponding to a subsequence of a read being mapped, that is built upby adding nucleotides to the subsequence, using the added nucleotides toselect next links to traverse in the tree, and going as deep asnecessary until a unique sequence has been generated. This uniquesequence may also be termed a seed, and may represent a branch and/orroot of the sequence tree data structure. Alternatively, the treedescent may be terminated before the accumulated subsequence is fullyunique, so that a seed may map to multiple locations in the referencegenome. Particularly, the tree may be built out for every startingposition for the reference genome, then the generated reads may becompared against the branches and/or roots of the tree and thesesequences may be walked through the tree to find where in the referencegenome the read fits. More particularly, the reads of the FASTQ file maybe compared to the branches and roots of the reference tree and oncematched therewith the location of the reads in the reference genome maybe determined. For example, a sample read may be walked along the treeuntil a position is reached whereby it is determined that theaccumulated subsequence is unique enough so as to identify that the readreally does align to a particular position in the reference, such aswalking through the tree until a leaf node is reached.

A disadvantage, however, of such a prefix and/or suffix tree is that itis a huge data structure that must be accessed a multiplicity of timesas the tree is walked so as to map the reads to the reference genome. Anadvantage of a hash table function, on the other hand, as described ingreater detail herein below, is that once built, it typically only takesone look up to determine where, if anywhere, there may be a matchbetween a seed and the reference. A prefix and/or suffix tree willtypically take a plurality of look ups, e.g., 5, 10, 15, 20, 25, 50,100, 1,000, or more, etc., in determining if and where there is a match.Further, due to the double helix structure of DNA, a reverse complementtree may also need to be built and searched, as the reverse complementto the reference genome may also need to be found. With respect to theabove, the data tree is described as being built from the referencegenome which is then compared with the reads from the subject'ssequenced DNA, however, it is to be understood that the data tree mayinitially be built from either the reference sequence or the samplereads, or both, and compared one to the other as described above.

Alternatively, or in addition to employing a prefix or a suffix tree, aBurrows/Wheeler transform can be performed on the data. For instance, aBurrows/Wheeler transform may be used to store a tree-like datastructure abstractly equivalent to a prefix and/or suffix tree, in acompact format, such as in the space allocated for storing the referencegenome. In various instances, the data stored is not in a tree-likestructure, but rather the reference sequence data is in a linear listthat may have been scrambled into a different order so as to transformit in a very particular way such that the accompanying algorithm allowsthe reference to be searched with reference to the sample reads so as toeffectively walk the “tree”. An advantage of the Burrows/Wheelertransform, such as over a prefix and/or suffix tree, is that ittypically requires less memory to store, and an advantage over a hashfunction is that it supports a variable seed length, and hence it can besearched until a unique sequence is determined and a match found. Forinstance, as with the prefix/suffix tree, however many nucleotides ittakes for a given sequence to be unique, or to map to a sufficientlysmall number of reference positions, determines the length of the seed.Whereas for a hash table, the seeds are all of the same predeterminedlength. A disadvantage, however, for the Burrows/Wheeler transform isthat it typically requires a multiplicity of lookups, such as two ormore look ups, such as for every step down the tree.

Alternatively, or in addition to utilizing one or both a prefix/suffixtree and/or a Burrows/Wheeler transform on the reference genome andsubject sequence data, so as to find where the one maps against theother, another such method involves the production of a hash table indexand/or the performance of a hash function. The hash table index may be alarge reference structure that is built up from sequences of thereference genome that may then be compared to one or more portions ofthe read to determine where the one may match to the other. Likewise,the hash table index may be built up from portions of the read that maythen be compared to one or more sequences of the reference genome andthereby used to determine where the one may match to the other.

More particularly, in any of the mapping algorithms described herein,such as for implementation in any of the method steps herein disclosed,one or all three mapping algorithms, or others known in the art, may beemployed, in software or hardware, so as to map one or more sequences ofa sample of sequenced DNA with one or more sequences of one or morereference genomes. As described herein in greater detail below, all ofthese operations may be performed via software or by being hardwired,such as into an integrated circuit, such as on a chip, for instance aspart of a circuit board. For instance, the functioning of one or more ofthese algorithms may be embedded onto a chip, such as into a FPGA (fieldprogrammable gate array) or ASIC (application specific integratedcircuit) chip, and may be optimized so as to perform more efficientlybecause of their implementation in such hardware.

Additionally, one or more, e.g., two or all three, of these mappingfunctions may form a module, such as a mapping module, that may formpart of a system, e.g., a pipeline, that is used in a process fordetermining an actual entire genomic sequence, or a portion thereof, ofan individual. The output returned from the performance of a mappingfunction may be a list of possibilities as to where one or more, e.g.,each, read maps to one or more reference genomes. For instance, theoutput for each mapped read may be a list of possible locations the readmay be mapped to a matching sequence in the reference genome. In variousembodiments, an exact match to the reference for at least a piece, e.g.,a seed of the read, if not all of the read may be sought. Accordingly,in various instances, it is not necessary for all portions of all thereads to match exactly to all the portions of the reference genome.

Further, one or all of these functions may be programmed in such amanner that exact or approximate matching and/or editing, such asediting of the results, may be performed. Hence, all of these processescan be configured to do inexact matching as well, where desired, such asin accordance with a preselected variance, such as 80% matching, 85%matching, 90% matching, 95% matching, 99% matching, or more. However, asdescribed in greater detail herein below, inexact matching may be a lotmore expensive such as in time and processing power requirements,because it may require any number of edits, e.g., where the edit may bea SNP or insertion or deletion of one or more bases, e.g., 1 or 2 or 3or 5 or more edits, to be performed so as to achieve an acceptablematch. Such editing is likely to be used more extensively inimplementing hashing protocols or when implementing prefix and/or suffixtrees and/or performing a Burrows/Wheeler transform.

With respect to hash tables, a hash table may be produced in manydifferent ways. In one instance, a hash table may be built by breakingthe reference genome into segments of standard length, e.g., seeds ofabout 16 to about 30 nucleotides or more in length, such as about 18 toabout 28 nucleotides, formatting them into a searchable table, andmaking an index of all the reference segments from which sequenced DNA,e.g., one or more reads, or a portion thereof, may be compared todetermine matching. More particularly, a hash table index may begenerated by breaking down the reference genome into segments ofnucleotide sequences of known, uniform length, e.g., seeds, and storingthem in random order into individual cubicles in the reference table.This may be done for a portion or the entire reference genome so as tobuild an actual reference index table that may be used to compareportions of the reference genome with portions of one or more reads,such as from a FASTQ file, for the purpose of determining matching.

This method may then be repeated in approximately the same manner for aportion, e.g., a majority or all, of the reads in the FASTQ file, so asto generate seeds of the appropriate, e.g., selected, length. Forinstance, the reads of the FASTQ file may be used to produce seeds of apredetermined length, which seeds may be converted into binary form andfed through a hash function and fit into a hash table index where thebinary form of the seeds may match up with the binary segments of thereference genome, so as to give the location as to where in the genomethe sample seeds match with the position in the reference genome.

For example, where the read is approximately 100 bases long, a typicalseed may be about half or a about a third, e.g., about 27 to about 30bases, as long. Hence, in such an instance, for each read a multiplicityof seeds, e.g., approximately 3 or 4 seeds dependent on the length ofthe read and/or the length of the seeds, may be generated to cover theread. Each seed may then be converted into a binary form and/or then befed into the hash table and a possible result as to its position withrespect to the reference genome may be obtained. In such instances, theentire read need not be compared to every possible position in theentire reference genome, rather only a portion of the reads, e.g., oneor more of the generated sample seeds per read, need only be comparedsuch as to an index containing equivalent seed portions of the referencegenome. Hence, in various instances, a hash table may be configured suchthat by only one memory look up it can typically be determined where thesample seed and therefore read is positioned relative to the referencegenome. However, in certain instances, it may be desirable to perform ahash function and look up on one or more overlapping sections of seedsfrom one read. In such instances, the seeds to be generated may beformed in such a manner that at least a portion of their sequenceoverlaps one another. This may be useful for instance in getting aroundmachine and/or human errors or differences between the subject and thereference genome and may promote exact matching.

In certain instances, the building of the hash table as well as theperformance of one or more of the various comparisons is executed by thehash function. The hash function is in part a scrambler. It takes aninput and gives what appears to be a random order to it. In thisinstance, the hash function scrambler breaks down the reference genomeinto segments of a preselected length and places them randomly in thehash table. The data may then be stored evenly across the whole storagespace. Alternatively, the storage space may be segmented and/or storagetherein may be weighted differently. More particularly, the hashfunction is a function that takes any input and gives a number, such asa binary pattern out, which number may typically random except that forany one given input the same output is always returned. Hence, even iftwo inputs that are fed into the hash table are almost the same, becausethey are not an exact match, two completely, randomly different outputswill be returned.

Further, since genetic material may be composed of four basicnucleotides, e.g., “A”, “C”, “G”, and “T” (or “U” in the case of RNA),the individual nucleotides of the sequences, e.g., the referencesegments and or reads, or portions thereof, to be fed into the hashtable may be digitized and represented in binary format, such as whereeach of the four bases represents a two bit digital code, e.g., “A”=00,“C”=01, “G”=11, and “T”/“U”=10. In certain instances, it is this binary“seed” value that is then randomly placed in the hash table at a knownlocation having a value equal to its binary representation. The hashfunction, therefore, works to break down the reference genome intobinary representations of reference seeds and inserts each binary seeddata into a random space, e.g., cubicle, in the hash table based on itsnumeric value. Along with this digital binary code, e.g., access key,each cubicle may also include the actual entry points to where thesegment originated from in the actual reference genome, e.g., thereference position. The reference position therefore may be a numberindicating the position of the original reference seed in the genome.This may also be done for overlapping positions, which are put into thetable in random order but at known location, such as by the hashfunction. In a manner such as this, a hash table index may be generated,wherein the index includes the digital binary code for a portion or allof a plurality of segments of one or more reference genomes, which maythen be referenced by one or more sequences of genetic material, e.g.,one or more reads, or portions thereof, from one or more individuals.

When implementing the hash table and/or function as a module, such as amodule in a pipeline of modules, on software (such as where the bitwidth is 2× the number of bases in the seed described above) and/orhardware, as referenced above, the hash table can be built so that thebinary representation of the reference seeds can be any bit widthdesired. As the seeds can be long or short, the binary representationscan be greater or lesser, but typically the seed length should be chosenso as to be long enough to be unique, but not too long that it is toohard to find matches between the seeds of the genome reference and theseeds of the sample reads, such as because of errors or variants. Forinstance, as indicated above, the human genome is made up of about 3.1billion base pairs, and a typical read may be about 100 nucleotides inlength. Hence, a useful seed length may be between about 16 or about 18nucleotides or less in length to about 28 or about 30 nucleotides ormore in length. For example, in certain instances, the seed length maybe a segment of 20 nucleotides in length. In other instances, the seedlength may be a segment of 28 nucleotides in length.

Consequently, where the seed length is a segment of 20 nucleotides, eachsegment may be represented digitally by a 40 bit output, e.g., a 40 bitbinary representation of the seed. For example, where 2 bits areselected to represent each nucleotide, e.g., such as where A=00, C=01,G=10, and T=11, a seed of 20 nucleotides×2 bits per nucleotide=a 40 bit(5 byte) vector, e.g., number. Where the seed length may be 28nucleotides in length, the digital, e.g., binary, representation of theseed may be a 56 bit vector. Hence, where the seed length isapproximately 28 nucleotides in length, 56 bits can be employed tohandle a 28 nucleotide seed length. More particularly, where the 56 bitsrepresents the binary form of the seeds of the reference genome thathave been randomly positioned in the hash table, a further 56 bits canbe used to digitally represent the seeds of the read that are to bematched against the seeds of the reference. These 56 bits may be runthrough a polynomial that converts the 56 bits in to 56 bits out in a1:1 correspondence. Without increasing or decreasing the number of bitsof output, performing this operation randomizes the storage location ofadjacent input values so that the various seed values will be uniformlydistributed among all possible storage locations. This also serves tominimize collisions among values that hash to the same location. Inparticular, in a typical hash table implementation described herein,only a portion of the 56 bits is used as a lookup address to select astorage location and the remaining bits are stored in that location forconfirmation of a match. If a hashing function were not used, a greatmany patterns having the same address bits, but different stored bitswould have to share the same hash location.

More specifically, there is similarity between the way the hash table isconstructed, e.g., by software and/or hardware placing the referencegenome seeds randomly in the hash table, and the way the hash table isaccessed by the seeds of the reads being hashed such that they bothaccess the table in the same way. Hence, seeds of the reference andseeds of the sample read that are the same, e.g., have the same binarycode, will end up in the same location, e.g., address, in the tablebecause they access the hash table in the same manner, e.g., for thesame input pattern. This is the fastest known method for performing apattern match. Each lookup takes a nearly constant amount of time toperform. This may be contrasted with a Burrows-Wheeler method which mayrequire many probes (the number may vary depending on how many bits arerequired to find a unique pattern) per query to find a match, or abinary search method that takes log₂(N) probes where N is the number ofseed patterns in the table.

Further, even though the hash function can break the reference genomedown into segments of seeds of any given length, e.g., 28 base pairs,and can then convert the seeds into a digital, e.g., binary,representation of 56 bits, not all 56 bits need be accessed entirely atthe same time or in the same way. For instance, the hash function can beimplemented in such a manner that the address for each seed isdesignated by a number less than 56 bits, such as about 20 to about 45bits, such as about 25 to about 40 bits, such as about 28 to about 35bits, including about 28 to about 30 bits may be used as an initial keyor address so as to access the hash table.

For example, in certain instances, about 26 to about 29 bits may be usedas a primary access key for the hash table, leaving about 27 to about 30bits left over, which may be employed as a means for double checking thefirst key, e.g., if both the first and second keys arrive at the samecell in the hash table, then it is relatively clear that said locationis where they belong. Specifically, in order to save space and reducethe memory requirements and/or processing time of the hash module, suchas when the hash table and/or hash function are implemented in hardware,the about 26 to about 29 bits representing the primary access keyderived from the original 56 bits representing the digitized seed of aparticular sequenced read may be employed by the hashing function tocomprise the primary address, leaving about 27 to about 30 bits that canbe used in a double checking method.

More particularly, in various instances, about 26 to about 29 bits fromthe 56 bits representing the binary form of a reference seed may beemployed to comprise a primary address, which designated 26 to 29 bitsmay then be given a randomized location in the hash table, which in turnmay then be populated with the location of where the reference seedoriginally came from along with the remaining 27 to 30 bits of the seedso that an exact match may be ascertained. The query seeds representingthe reads of the subject genome converted into binary form may also behashed by the same function in such a manner that they as well arerepresented by 29 bits comprising a primary access key. If the 29 bitsrepresenting the reference seed are an exact match to the 29 bitsrepresenting the query seeds, they both will be directed to the sameposition in the hash table. If there was an exact match to the referenceseed, then we expect to find an entry at that location containing thesame remaining 27 to 30 bits. In such an instance, the 29 designatedaddress bits of the reference sequence may then be looked up to identifythe position in the reference to where the query read from which thequery seed was derived, aligns.

However, with respect to the left over 27 to 30 bits, these bits mayrepresent a secondary access key that may also be imported into the hashtable as well, such as for the purpose of ensuring the results of thefirst 26 to 29 bits of the primary access key. Because the hash tablerepresents a perfect 1:1 scrambling of the 28 nucleotide/56 bitsequence, and only about 26 to about 29 of the bits are used todetermine the address, these 26 to 29 bits of the primary access keyhave basically been checked, thereby determining the correct address ina first go around. This data, therefore, does not need to be confirmed.However, the remaining about 27 to about 30 bits of the secondary accesskey must be checked. Accordingly, the remaining about 27 to 30 bits ofthe query seeds are inserted into the hash table as a means forcompleting the match. Such an implementation may be shorter than storingthe 56 bit whole key, and thus, saves space and reduces over all memoryrequirements and processing time of the module.

The hash table, therefore, can be configured as an index where knownsequences of one or more reference genomes that have been broken downinto sequences of predetermined lengths, e.g., seeds, such as of 28nucleotides in length, are organized into a table randomly, and one ormore sequenced reads, or “seed” portions thereof, derived from thesequencing of a subject's genomic DNA or RNA, may be passed through thehash table index, such as in accordance with a hash function, so as tolook up the seed in the index, and one or more positions, e.g.,locations in the reference genome, may be obtained from the table wherethe sample seed matches positions in the reference genome. Using a bruteforce linear search to scan the reference genome for locations where aseed matches, over 3 billion locations would have to be checked.However, by using a hashing approach, each seed lookup can occur inapproximately a constant amount of time. Often, the location can beascertained in a single access. In cases where multiple seeds map to thesame location in the table, a few additional accesses may be made tofind the seed being currently looked up. Hence, even though there can be30M or more possible locations for a given 100 nucleotide length read tomatch up to, with respect to a reference genome, the hash table and hashfunction can quickly determine where that read is going to show up inthe reference genome. By using a hash table index, therefore, it is notnecessary to search the whole reference genome to determine where theread aligns.

As indicted above, chromosomes have a double helix structure that iscomprised of two opposed, complementary strands of nucleic acidsequences that are bound together so as to form the double helix. Forinstance, when the double helix structure is formed these complementarybase pairs bind one with the other in accordance with the followingformula: “A” binds to “T”, and “G” binds to “C”. Accordingly, thisresults in two equal and opposite strands of nucleic acid sequences thatare the complement of each other. More particularly, the bases of anucleotide sequence of one strand will be mirrored by theircomplementary bases on the opposed strand resulting in two complementarystrands. However, transcription of DNA takes place in one directiononly, starting from one end of the DNA and moving towards the other.Hence, as it turns out, for one strand of the DNA, transcription takesplace in one direction, and for its complement strand, transcriptiontakes place in the opposite direction. Consequently, the two strands ofDNA sequences turn out to be reverse complemented, that is if thesequence order of one strand of the DNA is compared to the other whatcan be seen is two strands where the nucleotide letters of one strandare switched for their complement in the other strand, e.g., “As” for“Ts” and “Gs” for “Cs” and vice versa, and their order is reversed.

Because of the double helix structure of the DNA, during the sample prepstep prior to sequencing the DNA, the chromosomes are pulled apart,e.g., de natured, separated into separate strands, and then lysed intosmaller segments of a predetermined length, e.g., of 100-300 bases long,which are then sequenced. It is possible to separate the strands priorto sequencing so that only one strand is sequenced, but typically thestrands of DNA are not separated and so both strands of DNA aresequenced. Accordingly, in such an instance, about half of the reads inthe FASTQ file may be reverse complemented.

Of course, both strands of the reference genome, e.g., the complementand the reverse complement, may be processed and hashed as describedabove, however this would make the hash table twice as big, and make theperformance of the hash function take twice as long, e.g., it couldrequire about twice the amount of processing to compare both complementand reverse complemented sequences of the two genomic sequences.Accordingly, to save memory space, reduce processing power, and/ordecrease the time of processing, in various instances, only one strandof the model genomic DNA need be stored in the hash table as areference.

However, because in accordance with typical sequencing protocols, suchas where the two strands of the subject DNA have not been isolated fromone another, any read generated from the sequenced DNA can be fromeither strand, the complement or its reverse complement, it may bedifficult to determine which strand is being processed, the complementof the reverse complement. More specifically, in various instances,since only one strand of the reference genome need be used to generatethe hash table, half of the reads generated by the sequencing protocolmay not match the particular strand, e.g., either the complement or itsreverse complement, of the model genome reference, e.g., because halfthe time the read being processed is a reverse complement with respectto the hashed segments of the reference genome. Hence, only the readsgenerated from one strand of the DNA will match the indexed sequences ofthe reference genome, while the reads generated from the other strandwill theoretically be their reverse complements and will not matchanywhere in the reference genome. Further, an additional complicationcan be that for any given read that is reverse complemented to thestored reference genome strand, the read may still, erroneously, matchto a portion of the reference genome, such as by mere chance. In view ofthe above, in order for mapping to proceed efficiently, in variousinstances, it not only must be determined where the read matches in thereference genome it must also be determined if the read is reversecomplemented. Therefore, the hash table and/or function module should beconstructed so as to be able to minimize these complications and/or thetypes of errors that may result therefrom.

For instance, as indicated above, in one instance, the hash table couldbe populated with both the complement and the reverse complement for thereference genome so that every read or its reverse complement of thesubject's sequenced DNA can be matched to its respective strand in thegenomic reference DNA. In such an instance, for any given seed in aread, the seed should theoretically match with one strand or the other,the complement or the reverse complement of the reference, assuming noerrors or variations. However, storing both strands of the referencegenome in the hash index can require about twice as much storage space(e.g., instead of 32 gigabytes 64 gigabytes may be necessary), and mayrequire twice the amount of processing resources and/or twice as muchtime for processing. Further, such a solution doesn't solve the problemof palindromes that can match in both directions, e.g., the complementand reverse complement strands.

Accordingly, although the hash table index may be constructed to includeboth strands of the genomic reference sequence. In various instances,the hash table may be constructed so as to only include one strand ofthe model genome as a reference. This may be useful because storing thehash table in memory will require half of the storage and/or processingresources than would be required if both strands were to be stored andprocessed, and thus, the time required for a look up should also requireless time. However, storing only one strand of the genome as a referencecould cause complications because, as indicated above, where thesequenced subject DNA is double stranded, it is not typically known fromwhich strand any given read was generated. In such an instance,therefore, the hash table should be constructed to account for the factthe read being mapped may be from either strand and thus can be thecomplement or reverse complement of the stored segments of the referencegenome.

Accordingly, in various instances, such as where only one orientation ofseeds from the reference are populated into the hash table, whenperforming the hash function on the seeds generated from the reads ofthe FASTQ file, the seed may first be looked up in its presentorientation, and/or may then be reverse complemented and the reversecomplement may be looked up. This may require two looks up in the hashindex, e.g., twice as many, but one of the seed or its reversecomplement should match its complementary segment in the referencegenome, assuming no errors or variations, and it should reduce theoverall processing resources, e.g., less memory is used, as well asreducing time, e.g., not as many sequences need to be compared.

More particularly, such as where a seed in one particular orientation iscomprised of 28 nucleotides, e.g., digitally represented in a 56 bitbinary format, as described above, the seed can be reverse complementedand the reverse complement can also be represented digitally in a 56 bitbinary format. The binary format for each representation of the seedsequence and its complement results in a number, e.g., an integer,having a value represented by that number. These two values, e.g., thetwo integers, may be compared and the number with the higher or lowervalue, e.g., higher or lower absolute value, may be selected as thecanonical choice of orientation and that is the one that can be storedin the hash table and/or subjected to the hash function. For instance,in certain instances, the number with the higher value may be selectedfor being processed by the hash function.

Another method that may be employed is to construct seeds wherein eachseed is comprised of an odd number of bases. The canonical orientationto be selected then may be those strands having a middle base being an“A” or a “G”, but not a “T” or a “C”, or vice versa. The hash functionthen will be performed on the seeds meeting the requirements of thecanonical orientation. In such a manner, it is only the two bitsrepresenting the middle base that needs to be compared to see which hasthe higher value and it is only the 2 bits of that sequence that arelooked up. Hence, you only have to look at the bits representing themiddle two bases. Typically, this can work well because if the seed isan odd length, then it always reverse complements the center base.However, although this may work for odd seed lengths, hashing thoseseeds having a higher, or lower, value, as described above, should workfor all seed lengths, albeit such a method may require having toprocess, e.g., look up, more bits of data.

These methods may be performed for any number of seeds, e.g., all seedsof the reference and/or any number of seeds, e.g., all, derived from allor a portion of the reads of the FASTQ file. Approximately half of thetime the binary representation of the seeds of a given orientation,e.g., the complement, will have a higher value, and approximately halfthe time the binary representation of the seeds of the oppositeorientation, e.g., the reverse complement, will have the higher value.But, when looking at the binary numbers, whichever one has the highervalue, that is the one that gets fed into the hash table. For instance,the binary integers for each read and its complement may be compared,and the sequence having the first 1 encountered is the one of the twostrands selected to be stored as the strand in the hash table and/or besubjected to the hash function. If both strands have a first 1 in thesame position, then the strand having the second 1 that comes first isselected, and so on. Of course, the read with the lower value may alsobe selected, in which case the strand having the first and/or largernumber of initial 0's will be selected. An indication, e.g., a flag, mayalso be inserted into the hash table where the flag indicates whichorientation, complement or reverse complement, the stored and/or hashedstrand represents, e.g., a 1RC flag, if reverse complemented.

More particularly, when performing the hash function and accessing thehash table, seeds from the genomic reference DNA and seeds derived fromthe reads of the sequence data are subjected to these same operations,such as converted into binary form and compared with its reversecomplement where the integers having the higher, or lower, values areselected as the canonical orientations and subjected to the hashfunction and fed into the hash table to be looked up and matched againsteach other. However, because it is the same operation being performed insubstantially the same manner on the reference sequences and the readsequences, the same record will be derived, if the two sequences, thereference and the subject seeds, have the same sequence to begin with,even if one was reverse complemented, they will all be directed to thesame cell in the hash table.

Consequently, if a certain seed in the reference having a given sequencein a particular orientation is converted to binary form and hashed, andthen a seed derived from a sample read having the same sequence, but inits reverse orientation, e.g., reverse complemented, and it is subjectedto the above protocols, because of the above disclosed methods, when thebinary value is determined and the hash function performed, the look upwill be directed to the very same address in the hash table as if thehash function were performed on the complimentary seed to begin with.Hence, in this manner it doesn't matter which orientation the seed beingprocessed is in because it will always be directed to the same address.

Therefore, in a manner such as this, the methods herein disclosed areable to hash and thereby determine the location of the seed within thetable despite its orientation, and because of the flag in the record itwill also be known if any given seeds is reverse complemented. Forinstance, it will be known if the seed was flipped from the referenceand it will also be known if the seed derived from the subject read hadto be flipped as well. Consequently, if the decision was the same onboth sides then the orientation is the same between the read and thereference. However, if one side is flipped and the other is not, then itcan be concluded that the read maps reverse complemented to thereference. Hence, by using a hash table it may be determined where inthe genome a given read, or portion thereof, e.g., a seed, matchesand/or if it is reverse complimented. Further, it is to be understoodthat although the above is described with respect to generating the hashtable from the reference genome and performing various ancillary hashfunction processes on the seeds generated from the reads, e.g., from aFASTQ file, the system can also be structured such that the hash tableindex is generated from seeds derived from the reads of the subject'ssequenced DNA, and the various ancillary hash function processes, asherein described, are performed on seeds generated from the referencegenome.

As set forth above, an advantage of employing a hash table and/or a hashfunction is that by employing the use of seeds, a majority of the readsof the sequenced DNA can be matched to the reference genome often byemploying single hash lookups, and in various instances, not all seedsderived from a read need be hashed and/or looked up. Seeds may be of anysuitable length, such as relatively short, e.g., 16 nucleotides or less,such as about 20 nucleotides, such as about 24 nucleotides, such asabout 28 nucleotides, such as about 30 or about 40 or about 50, or 75 orabout 100 nucleotides, or even up to 250 or 500, or 750, or even 999 oreven about 1,000 nucleotides in length; or relatively long such as overabout 1,000 nucleotides or over about 10,000, or over about 100,000 orover 1,000,000 or more nucleotides in length. However, as describedabove, there are some disadvantages to using seeds, such as in a hashtable, in particular with respect to selecting seeds of the appropriatelength.

For instance, any suitable seed length may be employed in a mappingfunction, however there are advantages and disadvantages of usingrelatively short or relatively long seed lengths. For example, theshorter the seed length the less likely it is to incorporate an error ora variation that can prevent finding a match within the hash table.However, the shorter the seed length, the less unique it is, and themore matching is to be expected between the seeds of the referencegenome and the seeds derived from the reads of the subject's sequencedDNA. Further, the shorter the seed length the more lookups will have tobe performed by the hash function, taking more time and increasedprocessing power.

On the other hand, the longer the seed length the more unique it is andthe less likely there is to be multiple matching positions between theseeds between the seeds of the reference and the query. Also, with alonger seed, there need be fewer seeds within the read, so fewer lookups, thereby taking less time and requiring less processing power. Thelonger the seed, however, the more likely it is that the seeds derivedfrom the sequenced DNA may include an error, such as a sequencing errorand/or may incorporate a variation as compared to the reference thuspreventing a match from being made. Longer seeds further have thedisadvantage of being more likely to hit the end of the read and/or theend of the chromosome. Hence, where a seed is only 20-100 nucleotides inlength, there may be several matches within the hash table, however,where the seed is 1,000 or more nucleotides in length there may be muchfewer matches, but there may be no matches at all.

There are some methods for helping to minimize these issues. One methodis to ensure there is appropriate oversampling generated in the DNAprocessing steps prior to sequencing. For instance, as it is known thatthere is typically at least one variation within every 1,000 base pairs,the seed length may be chosen to maximize matches, while at the sametime minimizing non-matches due to the incorporation of errors and/orvariants. Additionally, the use of oversampling, such as in thepre-sequencing and/or sequencing steps, can be employed as a furthermethod for minimizing various problems that are inherent to using seeds,such as within a hash function.

As indicated above, oversampling produces pileups. Pileups are thosecollections of reads that map in an overlapping fashion generally to thesame place in the genome. For the majority of sample reads, such pileupsmay not be necessary, such as where the reads, and/or seeds generatedtherefrom, do not include a variant and/or do not map to multiplepositions in the hash table (e.g., are not exactly duplicated in thegenome). However, for those reads and/or seeds that may include avariant and/or an error and/or other mismatch between the seed and/orread and the reference genome, the production of pileups for any givenregion of the genome may be useful. For instance, even though only oneexact hit between a seed generated from a read of the sample genome isnecessary so as to be able to map the sample read to the referencegenome, however, the fact that there may be a machine error or a truevariant in the sample DNA sequence that could prevent such an exactmatch between the read and the reference from occurring, often timesmakes the production of overlapping pileups in the pre-sequencing andsequencing steps useful.

For example, for those instances where a sample seed does in factcontain a variant or an error, the production of read pileups may beuseful in distinguishing between actual variance and machine and/orchemistry errors. In such an instance, a pileup can be employed todetermine whether an apparent variation is in fact a real variation. Forinstance, if 95% of the reads in the pileup indicate that there is a “C”in a certain position, then odds are that is the correct call, even ifthe reference genome has a “T” at that location. In such an instance,the mismatch may be due to a SNP, e.g., a substitution of a “C” for a“T” in that position in the genome, where the genetic code for theindividual actually varies from that of the reference. In such aninstance, the depth of the pileup may be employed so as to compare theoverlapping portions of the reads of the pileup at a position wherethere is variance, and based on the percentage of reads in the pileuphaving the variance, it can be determined whether the variance is infact due to an actual variation in the sample sequence. Accordingly, theactual sequence of the reads that best fits the genomic sequence, may inpart be determined based on what is reflected in the pileup depths. Thedisadvantage of using pileups, however, is that it requires moreprocessing time to process all the excess reads and/or seeds generatedthereby.

Another method for minimizing the issues inherent in short or long readsis to employ a secondary hash table along with or in conjunction withthe first, e.g., primary hash table. For instance, a second hash tableand/or hash function may be employed for those seeds that do not haveany hits in the primary hash table, or for those seeds that havemultiple hits in the primary hash table. For example, when comparing oneseed with another there are several outcomes that may result. In oneinstance, a no hit, e.g., a no match anywhere between the two sequences,may result, in which case this suggests a possible error or variationsuch as in the seed of a read of the subject as compared against a seedderived from the reference genome. Or there may be one or a plurality ofmatches found. If a large number of matches are found, however, thiscould be problematic.

For instance, with respect to the primary hash table, if each seed inthe reference being hashed appears only a few times, e.g., once, twice,or three times, etc. then there may not be a need for a secondary hashtable and/or hash function. However, if one or more of the seeds occursa greater number of times, e.g., 5, 10, 15, 20, 25, 50, 100, 1,000, ormore times, this could be problematic. For example, there are knownregions in the sequence of the human genome that have been determined tobe mathematically significant in that they are repeated a multiplicityof times. Consequently, any seed mapping to one of these positions, mayin fact inadvertently map to a multiplicity of these positions, such aswhere the seed comprises the nucleotides of the overlapping sequences.In such an instance, determining which out of all the possibilities theseed actually aligns to may be difficult. However, as these repeatingregions are known, and/or become known, any seed that would typicallymap to one or more of these regions may be demarcated to be allocated toa secondary hash table for processing by the first or a secondary hashfunction, so as to not waste time and processing power trying to use aprimary hashing function to determine something that is likely to beindeterminable.

More particularly, when comparing the seeds of the genomic reference tothe seeds generated from the subject's genomic reads, anywhere from 1 tohundreds or even thousands of match positions may result. The presentsystem, however, may be configured to handle a certain number ofduplicative matches, such as without the need for further processingsteps, such as where the number of matches is below about 50, or belowabout 40, or below about 30, such as below about 25 or about 20, such asbelow about 16 matches or below about 10 or about 5 matches. However, ifthere are more matches of viable hits than this that are returned, thenthe system can be configured to implement a secondary hash function,e.g., using a secondary hash table.

Accordingly, rather than placing such seeds known to have an increasedlikelihood of redundancy in the primary hash table, such seeds can beplaced in a secondary hash table, or a secondary region in the firsthash table. Additionally, in some instances, a record that doesn'tcommunicate anything about the multiplicity of potential map positionsfor that seed, but rather communicates a command to access a secondaryhash table, e.g., an extend record, can be placed in the primary hashtable. For example, the extend record can be an instruction, such as aninstruction to extend the primary, e.g. non unique or duplicative, seedlength to a longer, more unique seed length, such as by adding on one ormore additional bases next to it, e.g., on the end(s) of the seed, tomake it a longer seed sequence that can then get hashed and looked up,such as in the secondary table.

The record can be configured such that it informs or otherwise instructshow much to extend the known redundant seed by a given amount, and mayalso instruct as to where and/or how to extend the seed. For instance,because the hash table is usually precomputed, e.g., originallyconstructed from the seeds generated from the reference genome(s), itmay be known prior to constructing the table, which, if any, of theseeds generated from the reference genome are going to occur amultiplicity of times. Hence, in various instances, it may bepredetermined which seeds are going to need to be shifted over to thesecondary hash table. For example, when constructing the hash tableindex, the characteristics of the reference seed sequences being inputinto the hash table as an index are known, so for every potential seedit may be determined whether it's a case that is going to give amultiplicity of hits, e.g., from 10-10,000 hits.

More particularly, in various instances, an algorithm can be performedto determine all the predicted matches a given seed derived from thereference and/or the subject's reads may have. If it is determined thatfor any particular seed that it is likely to return a multiplicity ofmatches, a flag, e.g., a record, may be generated, such as within a cellof the hash table, indicating that this particular seed is a highfrequency hit. In such an instance, the record can further instruct thatthe primary hashing of this seed, and such seeds like it, should beskipped over because it is not practical to perform the number, e.g.,20-10,000 or more evaluations on such a seed needed to accuratelydetermine where the seed actually maps. In such an instance, the primaryhash function may not be able to accurately determine which position outof all the possible positions to where the seed may match, is the one towhere the read actually aligns, and thus for practical purposes, becausethe seed cannot accurately be mapped at this stage, the primary hashfunction may not be likely to return a useable result, such as a resultindicating accurately where the seed actually matches in the genome.

In such an instance, the hash function algorithm may be configured tocalculate what would need to be done to make the redundant seed moreunique. For example, the secondary hash function may determine by howmany bases the seed needs to be extended, and in what order, and in whatlocation, so as to ensure that the seed is no longer redundant, butrather suitably unique so as to be hashed. Accordingly, the record mayalso include an instruction to extend the redundant seed, e.g., extendby two, by four, by six, etc., on one or both ends of the seed so as toachieve a predetermined level of uniqueness. In such a manner as this,seeds that at first appear to be identical can be determined to benon-identical.

For example, in some instances, a typical record can instruct that theduplicative seed be extended by up to X number of odd or even bases, butin some instances, extended by an even number of bases, such as fromabout 2 to 4 to about 8 to 16 to about 32 or about 64 or more bases,such as equally on each side. For instance, where the extension is to beby 64 bases, the record could instruct that 32 bases be added on eachside of the seed. The number of bases by which the seed is to beextended is configurable and may be any suitable number dependent on howthe system is constructed. In certain instances, the secondary hashfunction may be employed to determine by how many bases the seed shouldbe extended so as to get a more reasonable number of match results back.Therefore, the extension may be to the point of relative uniqueness,such as to where there is only 1, 2, 3, or even up to 16 or 25 or 50match positions where the pattern shows up. In various instances,extending the seed equally from both ends may be useful such as to avoidproblems with reverse reads, but in various instances the seed may beextended by the addition of one or more bases unequally to both sides.

More particularly, such as in one example, if the seed includes 28bases, and an extend record, such as an extend record positioned withina cell in the primary hash table, instructs the hash function to extendthe seed, such as by 64 bases, then the record may further direct thehash function as to how to extend the seed, such as by adding 32 baseson each side of the seed. However, the extension can take place at anysuitable position on the read and may be done in a symmetrical orasymmetrical fashion. In certain instances, the record may instruct thehash function to extend the seed symmetrically because in certaininstances such a symmetrical extension may work better, such as withreverse complements, discussed herein. In such an instance, the samenumber of bases will be added such as to the opposite sides of the seedwhen extending. Although in other instances extension may be performedby adding an even or an odd number of bases in a non-symmetrical format,and hence, it is not necessary to extend the seed by same number ofbases on each side. Typically, the primary hash table is configured suchthat it is not completely full. For example it is desirable to configureit not to exceed 80% or 90% of its capacity. This is to maintain highperformance of the lookup rate. When there are a high number ofcollisions in hashing seeds to the same location when constructing thetable, the storing mechanism will create a chain of references to otherlocations so that the lookup mechanism will be able to find the oneassigned to the overflowed seed. The denser the table, the higher thenumber of collisions and the longer the chains to be followed to findthe actual match.

In various instances, such as where the initial, redundant seed is 28bases long, and the record instructs for it to be extended, such as from18 to 32 to 64 bases, such as on each opposed side of the seed, thedigital representation of the seed may be about 64 bases×2 bits perbase=128 bits. Accordingly, dependent on how the mapping module is setup, this may be too big for the primary hash table to process. Hence, incertain instances, to deal with the need for such extensive processing,in certain embodiments, the secondary hashing module can be configuredto store the information associated with larger seeds. Since the numberof seeds requiring extension is a fraction of the total number of seeds,the secondary hash table may be smaller than the primary hash table.However, in other instances, such as to reduce the processingrequirements of the module, e.g., to save bits, the known redundantportion of the sequence, e.g., the primary sequence, may be replaced bya preselected variable such as of a predetermined sequence length. Insuch an instance, since the redundant sequence is already known andidentified, it does not need to be digitally represented in itsentirety. Rather, in various instances, all that is really needed to bedone is to substitute the known, redundant sequence with a knownvariable sequence, and all that really needs to be looked up are theextension portions, e.g., wings, that have been added to either side ofthe variable sequence, since those are the only portions of the initialsequence that are non-redundant and new. Hence, in certain instances,the primary sequence may be replaced by a shorter unique identifier code(such as a 24 bit proxy instead of 56 bit representation) and then theextension bases can be added to the proxy, such as a 36 bit extension(e.g., totaling 60 bits) that can then be put into the extend record inthe primary table. In a manner such as this, the disadvantages of havingtoo short and/or too long of reads can be minimized and the benefit ofhaving only one or a few look ups in the hash table can be maintained.

As indicated above, the implementation of the above described hashfunction may be executed in software of hardware. An advantage ofimplementing the hash module in hardware is that the processes may beaccelerated and therefore performed in a much faster manner. Forinstance, where software may include various instructions for performingone or more of these various functions, the implementation of suchinstructions often requires data and instructions to be stored and/orfetched and/or read and/or interpreted, such as prior to execution. Asindicated above, however, and described in greater detail herein below,a chip can be hardwired to perform these functions without having tofetch, interpret, and/or perform one or more of a sequence ofinstructions. Rather, the chip may be wired to perform such functionsdirectly. Accordingly, in various aspects, the disclosure is directed toa custom hardwired machine that may be configured such that portions orall of the above described hashing module may be implemented by one ormore network circuits, such as integrated circuits hardwired on a chip,such as an FPGA or ASIC.

For instance, in various instances, the hash table index may beconstructed and the hash function may be performed on a chip, and inother instances, the hash table index may be generated off of the chip,such as via software run by a host CPU, but once generated it is loadedonto and employed by the chip, such as in running the hash module. Incertain instances, the chip may include any suitable number ofgigabytes, such as 8 gigabytes, such as 16 gigabytes, such as 32gigabytes, such as 64 gigabytes, such as about 128 gigabytes. In variousinstances, the chip may be configurable such that the various processesof the hash module are performed employing only a portion or all thememory resources. For example, where a custom reference genome may bebuilt, a large portion of the memory may be dedicated to storing thehash reference index and/or for storing reads and/or for reserving spacefor other functional modules to use, such as where 16 gigabytes arededicated to storing the reads, 8 gigabytes may be dedicated to storingthe hash index and another 8 gigabytes may be dedicated to otherprocessing functions. In another example, where 32 gigabytes arededicated to storing reads, 26 gigabytes may be dedicated for storingthe primary hash table, 2.5 gigabytes may be dedicated for storing thesecondary table, and 1.5 gigabytes may be dedicated for the referencegenome.

In certain embodiments, the secondary hash table may be constructed soas to have a digital presence that is larger than the primary hashtable. For instance, in various instances, the primary hash table can beconfigured to store hash records of 8 bytes each with 8 records per hashbucket totaling 64 bytes per bucket, and the secondary hash table can beconfigured to store 16 hash records totaling 128 bytes per bucket. Foreach hash record containing overflow hash bits matching the same bits ofthe hash key a possible matching position in the reference genome isreported. For the primary hash table therefore, up to 8 positions may bereported. For the secondary hash table up to 16 positions may bereported.

Regardless of being implemented in hardware or software, in manyinstances, it may be useful to structure the hash table to avoidcollisions. For instance, there may be multiple seeds that, because ofvarious system artifacts will want to be inserted into the hash table atthe same place regardless of whether there is a match there or not. Suchinstances are termed collisions. Often times, collisions can be avoided,in part, by the way the hash table is structured. Accordingly, invarious instances the hash table may be structured so as to avoidcollisions, and therefore may be configured to include one or morevirtual hash buckets.

In various instances, the hash table can be structured such that it isrepresented in an 8 byte, 16 byte, 32 byte, 64 byte, 128 byte format, orthe like. But in various exemplary embodiments it may be useful torepresent the hash table in a 64 byte format. This may be useful, forinstance, where the hash function is to make use of accessing a memory,such as a DRAM, e.g., in a standard DIMM or SODIMM form factor, such aswhere the minimum burst size is typically 64 bytes. In such an instance,the design of the processor for accessing a given memory will be suchthat the number of bytes needed to form a bucket in the hash table isalso 64, and therefore a maximized efficiency may be realized. However,if the table were to be structured in a 32 byte format, this would beinefficient because about half the bytes delivered in a burst wouldcontain information not needed by the processor. That would cut theeffective byte delivery rate in half. Conversely, if the number of bytesused to form a bucket in the hash table is a multiple of the minimumburst size, e.g., 128, there is no performance penalty as long as theprocessor actually needs all of the information returned in a singleaccess. Therefore, in instances where the optimal burst size of thememory access is at a given size, e.g., 64 bytes, the hash table can bestructured so burst size of the memory is optimally exploited, such aswhere the bytes allocated for representing bins in the hash table andprocessed by the mapping function, e.g., 64 bytes, are coincident withthe burst size of the memory. Consequently, where the memory bandwidthis a constraint, the hash table can be structured so as to optimallyexploit such constraints.

Further, it is to be noted, that although a record may be crammed into 8bytes, the hash function can be constructed such that it is not the casethat 8 bytes from the table are read so as to process one record, asthis could be inefficient. Rather, all 8 records in a bucket can be readat once, or some sub-portion thereof. This may be useful in optimizingthe processing speed of the system as, given the architecture describedabove, it would cost the same time at the same speed to process all 8records as it would for simply processing 1 record. Accordingly, incertain instances, the mapping module may include a hash table thatitself may include one or more subsections, e.g., virtual sections orbuckets, wherein each bucket may have 1 or more slots, such as 8 slots,such that one or more different records can be inserted therein such asto manage collisions. However, in certain circumstances, one or more ofsuch buckets may fill up with records, so a means may be provided forstoring additional records in other buckets and recording information inthe original bucket indicating that the hash table lookup mechanismneeds to look further to find a match.

Hence, in certain instances it may also be useful to employ one or moreadditional methods such as for managing collisions, one such method mayinclude one or more of linear probing and/or hash chaining. Forinstance, if it is not known what exactly is being searched in the hashtable or a portion thereof, such as in one bucket of the hash table, andthe particular bucket is full, then the hash lookup function can beconfigured such that if one bucket is full and is searched and thedesired record not found, then the function can be directed to step tothe next bucket, e.g., the +1 bucket, and that bucket can then bechecked. In such a manner, all buckets can be searched when looking fora particular record. Such searching, therefore, can be performedsequentially looking through one bucket to another until what is beinglooked for is found or it becomes clear that it is not going to befound, such as where an empty slot in at least one of the buckets isfound. Particularly, where each bucket is filled sequentially, and eachbucket is searched according to the sequence of filling, if an emptyslot is found, such as when searching sequentially through bucketslooking for a particular record, then the empty slot could be indicativeof the record not existing, because if it did exist, it would at leasthave been positioned in the empty slot, if not in the preceding buckets.

More particularly, where 64 bytes are designated for storing theinformation in a hash bucket wherein 8 records are contained, uponreceiving a fetched bucket, the mapping processor can operate on all 8records simultaneously to determine which are matches and which are not.For instance, when performing a look up such as of a seed from a readobtained from the sequenced sample DNA against a seed generated from thereference genome, the digital representation of the sample seed can becompared against the reference seeds in all, e.g., 8, records so as tofind a match. In such an instance, several outcomes may result. A directmatch may be found. A sample seed may go into the hash table and, insome instances, no match is found, e.g., because it is just not exactlythe same as any corresponding seed in the reference, such as becausethere was a machine or sequencing error with respect to that seed or theread from which it is generated, or because the person has a geneticsequence that is different from the reference genome. Or a the seed maygo into the hash table and a plurality of matches may be returned, suchwhere the sample seed matches to 2, 3, 5, 10, 15, 20, or more places inthe table. In such an instance, multiple records may be returned allpointing to various different locations in the reference genome wherethat particular seed matches, the records for these matches may eitherbe in the same bucket, or a multiplicity of buckets may have to beprobed to return all of the significant, e.g., match, results.

In certain instances, such as where space may become a limiting factorin the hash table, e.g., in the hash table buckets, an additionalmechanism for resolving collisions and/or for saving space mayimplemented. For instance, when space becomes limited, such as when morethan 8 records need to be stored in a bucket, or when for otherinstances it is desirable, a hash chaining function may be performed.Hash chaining can involve, for example, replacing a record containing aspecific position location in the genomic sequence with a recordcontaining a chain pointer that instead of pointing to a location in thegenome points to some other address, e.g., a second bucket in thecurrent hash table e.g. a primary or a secondary hash table. This hasthe advantage over the linear probing method of enabling the hash lookupmechanism to directly access the bucket containing the desired recordrather than checking buckets sequentially in order.

Such a process may be useful given the system architecture. Forinstance, the primary seeds being hashed, such as in a primary lookup,are positioned at a given location in the table, e.g., their originalposition, whereas the seeds being chained are being put in a positionthat may be different from their original bucket. Hence, as indicatedabove, a first portion of the digitally represented seed, e.g., about 26to about 29 bits, can be hashed and may be looked up in a first step.And, in a second step, the remaining about 27 to about 30 bits can beinserted into the hash table, such as in a hash chain, as a means forconfirming the first pass. Accordingly, for any seed, its originaladdress bits may be hashed in a first step, and the secondary addressbits may be used in a second, confirmation step. Hence, the firstportion of the seeds can be inserted into primary record location, andthe second portion may be fit into the table in secondary record chainlocation. And, as indicated above, in various instances, these twodifferent record locations may be positionally separated, such as by achain format record. Therefore, in any destination bucket of chaining achain format record may positionally separate the entries/records thatare for local primary first bucket accesses and probing and thoserecords that are for the chain.

Such hash chains can be continued for a multiplicity of lengths. Anadvantage of such chaining is that where one or more of the bucketsinclude one or more, e.g., 2, 3, 4, 5, 6, or more empty record slots,these empty slots can be used to store the hash chain data. Accordingly,in certain instances, hash chaining may involve starting with an emptyslot in one bucket and chaining that slot to another slot in anotherbucket, where the two buckets may be at remote locations in the hashtable. Additional care may be taken to avoid confusion between recordsplaced in a remote bucket as part of a hash chain, and “native” recordsthat hash directly into the same bucket. As usual, the remaining about27 to about 30 bits of the secondary access key are checked againstcorresponding about 27 to 30 bits stored in the records placed remotelyin the chained bucket, but due to the distant placement of the chainedbucket from the original hash bucket, confirming these about 27 to 30bits would not be enough to guarantee that a matching hash recordcorresponds to the original seed reaching this bucket by chaining, asopposed to some other seed reaching the same bucket by direct access.(e.g., confirming the about 27 to 30 bits may be a full verificationwhen the about 26 to 29 bits used for hash table addressing areimplicitly checked by proximity to the initial hash bucket accessed.)

To prevent retrieving a wrong hash record without needing to storeentire hash keys in the records, a positional system may be used in achained bucket. Accordingly, a chained bucket must contain a chaincontinuation format record, which contains a further chain pointer tocontinue the bucket chain if required; this chain continuation recordmust appear in a slot of the bucket after all “native” recordscorresponding to direct hash access, and before all remote recordsbelonging to the chain. During queries, before following any chainpointer, any records appearing after a chain continuation record shouldbe ignored, and after following any chain pointer, any records appearingbefore a chain continuation record should be ignored.

For example, where the buckets are about 75%-85% full, 8 buckets may bescanned and only 15-25 slots may be found that can be used, whereas withhash chaining these slots may be found over 2 or 3 or 4 buckets. In suchan instance, the number of probe or chain steps required to store a hashrecord matters because it influences the speed of the system. At runtime, if probing is necessary to find the record, a multiplicity of hashlook up accesses, e.g., a 64 byte bucket read, may need to be performedwhich slows the system down. Hash chaining helps to minimize the averagenumber of accesses that have to be performed, because more excess hashrecords can generally be populated per chained bucket, which can beselected from a wide region, than per probing bucket, which must besequentially next. Therefore, a given number of excess hash records cantypically be populated into a shorter sequence of chained buckets thanthe necessary sequence of probing buckets, which likewise limits thenumber of accesses required to locate those excess records in a query.Nevertheless, probing remains valuable for smaller quantities of excesshash records, because probing does not require a bucket slot to besacrificed for a chain pointer.

For example, after it has been determined where all the possible matchesare for the seeds against the reference genome, it must be determinedwhich out of all the possible locations a given read may match to is infact the correct position to which it aligns. Hence, after mapping theremay be a multiplicity of positions that one or more reads appear tomatch in the reference genome. Consequently, there may be a plurality ofseeds that appear to be indicating the exact same thing, e.g., they maymatch to the exact same position on the reference, if you take intoaccount the position of the seed in the read.

The actual alignment, therefore, must be determined for each given read.This determination may be made in several different ways. In oneinstance, all the reads may be evaluated so as to determine theircorrect alignment with respect to the reference genome based on thepositions indicated by every seed from the read that returned positioninformation during the hash lookup process. However, in variousinstances, prior to performing an alignment, a seed chain filteringfunction may be performed on one or more of the seeds. For instance, incertain instances, the seeds associated with a given read that appear tomap to the same general place as against the reference genome may beaggregated into a single chain that references the same region. All ofthe seeds associated with one read may be grouped into one or more seedchains such that each seed is a member of only one chain. It is suchchain(s) that then cause the read to be aligned to each indicatedposition in the reference genome. Specifically, in various instances,all the seeds that have the same supporting evidence indicating thatthey all belong to the same general location(s) in the reference may begathered together to form one or more chains. The seeds that grouptogether, therefore, or at least appear as they are going to be near oneanother in the reference genome, e.g., within a certain band, will begrouped into a chain of seeds, and those that are outside of this bandwill be made into a different chain of seeds.

Once these various seeds have been aggregated into one or more variousseed chains, it may be determined which of the chains actuallyrepresents the correct chain to be aligned. This may be done, at leastin part, by use of a filtering algorithm that is a heuristic designed toeliminate weak seed chains which are highly unlikely to be the correctone. Generally, longer seed chains, in terms of length spanned withinthe read, are more likely to be correct, and furthermore, seed chainswith more contributing seeds are more likely to be correct. In oneexample, a heuristic may be applied wherein a relatively strong“superior” seed chain, e.g. long or having many seeds, filters out arelatively weak “inferior” seed chain, e.g. short or having few seeds.In one variation, the length of an inferior chain determines a thresholdlength, e.g. twice as long, such that a superior chain of at least thethreshold length can filter it out. In another variation, the seed countof an inferior chain determines a threshold seed count, e.g. five timesas many seeds, such that a superior chain of at least the threshold seedcount can filter it out. In another variation, the length of an inferiorchain determines a threshold seed count, e.g. two times the seed countminus the seed length, such that a superior chain of at least thethreshold seed count can filter it out. In some variations, such as whenchimeric alignments of reads are desired, only superior seed chainssubstantially overlapping inferior seed chains within the read mayfilter them out.

This process weeds out those seeds that have a low probability of havingidentified a region of the reference genome where a high qualityalignment of the read can be found. It, therefore, may be useful becauseit reduces the number of alignments that need to be performed for eachread thereby accelerating the processing speed and saving time.Accordingly, this process may be employed, in part, as a tuning feature,whereby when greater speed is desired, e.g., high speed mode, moredetailed seed chain filtering is performed, and where greater overallaccuracy is desired, e.g., enhanced accuracy mode, less seed chainfiltering is performed, e.g., all the seed chains are evaluated.

In various embodiments, seed editing may be performed, such as prior toa seed chain filtering step. For instance, for each read, if all of theseeds of that read are subjected to a mapping function and none of themreturned a hit, then there may be a high probability that there was oneor more errors in the read, for instance, an error that the sequencermade. In such an instance, an editing function, such as a one-changeediting process, e.g., an SNP editing process, can be performed on eachseed, such as where a no match outcome was returned. For example, atposition X, a one change edit function may instruct that the designatednucleotide be substituted for one of the other 3 nucleotides and it isdetermined whether a hit, e.g., a match, is obtained by making thatchange, e.g., a SNP substitution. This one-change editing may beperformed in the same manner on every position in the seed and/or onevery seed of the read, e.g., substituting each alternative base foreach position in the seed. Additionally, where one change is made in oneseed, the effects that change would have on every other overlapping seedmay be determined in view of that one change.

Such editing may also be performed for inserts, such as where one of thefour nucleotides is added at a given insert position, X, and it isdetermined if a hit was obtained by making the substitution. This may bedone for all four nucleotides and/or for all positions (X, X+1, X+2,X+3, etc.) in the seed and/or all the seeds in the reads. Such editingmay also be performed for deletions, such as where one of the fournucleotides is deleted at a given position, X, in the seed, and it isdetermined if a hit was obtained by making the deletion. This may thenbe repeated for all positions X+1, X+2, X+3, etc. Such editing, however,can result in a lot of extra processing work and time, such as byrequiring a multiplicity of additional lookups, such as 2, or 3, or 4,or 5, or 10, or 50, or 100, or 200, etc. Nevertheless, such extraprocessing and time may be useful if by such editing an actual hit canbe determined, e.g., a match made, where before there was no match. Insuch an instance, it can then typically be determined that an error wasmade and further that it was corrected, thereby salvaging the read.

Additionally, a further heuristic may be employed so as to determinewhether an editing function should be performed or not, whereby thealgorithm performs a calculation to determine the probability that a hitwill be obtained if such editing were to be performed. If a certainthreshold probability is met, such as 85% likelihood, then such seedchain editing may be performed. For instance, the system can generatevarious statistics on the seed chains, such as calculating how many highfrequency hits are present and/or how many seed chains contain highfrequency hits, and thereby determine if seed chain editing is likely tomake a difference in determining matches. For example, if it isdetermined that there are a large proportion of high frequency hits,then, in such an instance, seed chain editing may be skipped because itis unlikely to make various of the sequences unique enough to give a hitwithin a reasonable number of hash table look ups, such as 100 or fewer,50 or fewer, 40 or fewer, 30 or fewer, 20 or fewer, or 10 or fewer. Suchstatistics can be reviewed and it may then be determined whether to doseed editing or not. For instance, if the statistics show that for anyone read, if half the positions show no match, and the others show highfrequency matches, then it is probably worth doing seed editing, becausewhere no matches are returned, there is probably an error, but if a lotof high frequency matches are returned it may simply not be worthperforming seed editing.

The outcome from performing one or more of these mapping, filtering,and/or editing functions is a list of reads which includes for each reada list of all the possible locations to where the read may matchup withthe reference genome. Hence, a mapping function may be performed so asto quickly determine where the reads of the FASTQ file obtained from thesequencer map to the reference genome, e.g., to where in the wholegenome the various reads map. However, if there is an error in any ofthe reads or a genetic variation, you may not get an exact match to thereference and/or there may be several places one or more reads appear tomatch. It, therefore, must be determined where the various readsactually align with respect to the genome as a whole.

Accordingly, after mapping and/or filtering and/or editing, the locationpositions for a large number of reads have been determined, where forsome of the individual reads a multiplicity of location positions havebeen determined, and it now needs to be determined which out of all thepossible locations is in fact the true or most likely location to whichthe various reads align. Such aligning may be performed by one or morealgorithms, such as a dynamic programming algorithm that matches themapped reads to the reference genome and runs an alignment functionthereon.

An exemplary aligning function compares one or more, e.g., all of thereads, to the reference, such as by placing them in a graphical relationto one another, e.g., such as in a table, e.g., a virtual array ormatrix, where the sequence of one of the reference genome or the mappedreads is placed on one dimension or axis, e.g., the horizontal axis, andthe other is placed on the opposed dimensions or axis, such as thevertical axis. A conceptual scoring wave front is then passed over thearray so as to determine the alignment of the reads with the referencegenome, such as by computing alignment scores for each cell in thematrix.

The scoring wave front represents one or more, e.g., all, the cells ofthe matrix, or a portion of those cells, which may be scoredindependently and/or simultaneously according to the rules of dynamicprogramming applicable in the alignment algorithm, such asSmith-Waterman, and/or Needleman-Wunsch, and/or related algorithms. Forexample, taking the origin of the matrix (corresponding to the beginningof the read and/or the beginning of a reference window of the conceptualscoring wave front) to be at the top-left corner, first only thetop-left cell at coordinates (0,0) of the matrix may be scored, e.g., a1-cell wave front; next, the two cells to the right and below atcoordinates (0,1) and (1,0) may be scored, e.g., a 2-cell wave front;next the three cells at (0,2), (1,1), and (2,0) may be scored, e.g., a3-cell wave front. These exemplary wave fronts may then extenddiagonally in straight lines from bottom-left to top-right, and themotion of the wave front from step to step is diagonally from top-leftto bottom-right through the matrix. Alignment scores may be computedsequentially or in other orders, such as by computing all the scores inthe top row from left to right, followed by all the scores in the nextrow from left to right, etc. In this manner the diagonally sweepingdiagonal wave front represents an optimal sequence of batches of scorescomputed simultaneously or in parallel in a series of wave front steps.

For instance, in one embodiment, a window of the reference genomecontaining the segment to which a read was mapped is placed on thehorizontal axis, and the read is positioned on the vertical axis. In amanner such as this an array or matrix is generated, e.g., a virtualmatrix, whereby the nucleotide at each position in the read may becompared with the nucleotide at each position in the reference window.As the wave front passes over the array, all potential ways of aligningthe read to the reference window are considered, including if changes toone sequence would be required to make the read match the referencesequence, such as by changing one or more nucleotides of the read toother nucleotides, or inserting one or more new nucleotides into onesequence, or deleting one or more nucleotides from one sequence.

An alignment score, representing the extent of the changes that would berequired to be made to achieve an exact alignment, is generated, whereinthis score and/or other associated data may be stored in the given cellsof the array. Each cell of the array corresponds to the possibility thatthe nucleotide at its position on the read axis aligns to the nucleotideat its position on the reference axis, and the score generated for eachcell represents the partial alignment terminating with the cell'spositions in the read and the reference window. The highest scoregenerated in any cell represents the best overall alignment of the readto the reference window. In various instances, the alignment may beglobal, where the entire read must be aligned to some portion of thereference window, such as using a Needleman-Wunsch or similar algorithm;or in other instances, the alignment may be local, where only a portionof the read may be aligned to a portion of the reference window, such asby using a Smith-Waterman or similar algorithm.

The size of the reference window may be any suitable size. For instance,since a typical read may be from about 100 to about 1,000 nucleotideslong, the length of the reference window accordingly, in some instances,may be from about 100 to 1,000 nucleotides long or longer. However, insome instances, the length of the reads may be greater, and/or thelength of the reference window can be greater such as about 10,000,25,000, 50,000, 75,000, 100,000, 200,000 nucleotides long or more. Itmay be advantageous for the reference window to be padded somewhatlonger than the read, such as including 32 or 64 or 128 or 200 or even500 extra nucleotides in the reference window beyond the extremes of thereference genome segment to which the read was mapped, such as to permitinsertions and/or deletions near the ends of the read to be fullyevaluated. For instance, if only a portion of the read was mapped to asegment of the reference, extra padding may be applied to the referencewindow corresponding to the unmapped portions of the read, or longer bysome factor, such as 10% or 15% or 20% or 25% or even 50% or more, so asto allow the unmapped portions of the read space to fully align to thereference window. In some instances, however, the length of thereference window may be selected to be shorter than the length of thereads, such as where a long portion of the read is not mapped to thereference, such as more or less than 1000 nucleotides at one end of theread, such as in order to focus the alignment on the mapped portion.

The alignment wave front may be of unlimited length, or limited to anysuitable fixed length, or of variable length. For instance, all cellsalong the entire diagonal line of each wave front step extending fullyfrom one axis to the other axis may be scored. Alternatively, a limitedlength, such as 64 cells wide, may be scored on each wave front step,such as by tracing a diagonally 64-cell wide band of scored cellsthrough the matrix, and leaving cells outside of this band unscored. Insome instances, it may be unnecessary to calculate scores far from aband around the true alignment path, and substantial work may be savedby computing scores only in a limited bandwidth, using a fixed lengthscoring wave front, as herein described.

Accordingly, in various instances, an alignment function may beperformed, such as on the data obtained from the mapping module. Hence,in various instances, an alignment function may form a module, such asan alignment module, that may form part of a system, e.g., a pipeline,that is used, such as in addition with a mapping module, in a processfor determining the actual entire genomic sequence, or a portionthereof, of an individual. For instance, the output returned from theperformance of the mapping function, such as from a mapping module,e.g., the list of possibilities as to where one or more or all of thereads maps to one or more positions in one or more reference genomes,may be employed by the alignment function so as to determine the actualsequence alignment of the subject's sequenced DNA.

Such an alignment function may at times be useful because, as describedabove, often times, for a variety of different reasons, the sequencedreads do not always match exactly to the reference genome. For instance,there may be an SNP (single nucleotide polymorphism) in one or more ofthe reads, e.g., a substitution of one nucleotide for another at asingle position; there may be an “indel,” insertion or deletion of oneor more bases along one or more of the read sequences, which insertionor deletion is not present in the reference genome; and/or there may bea sequencing error (e.g., errors in sample prep and/or sequencer readand/or sequencer output, etc.) causing one or more of these apparentvariations. Accordingly, when a read varies from the reference, such asby an SNP or indel, this may be because the reference differs from thetrue DNA sequence sampled, or because the read differs from the true DNAsequence sampled. The problem is to figure out how to correctly alignthe reads to the reference genome given the fact that in all likelihoodthe two sequences are going to vary from one another in a multiplicityof different ways.

Accordingly, in various instances, the input into an alignment function,such as from a mapping function, such as a prefix/suffix tree, or aBurrows/Wheeler transform, or a hash table and/or hash function, may bea list of possibilities as to where one or more reads may match to oneor more positions of one or more reference sequences. For instance, forany given read, it may match any number of positions in the referencegenome, such as at 1 location or 16, or 32, or 64, or 100, or 500, or1,000 or more locations where a given read maps to in the genome.However, any individual read was derived, e.g., sequenced, from only onespecific portion of the genome. Hence, in order to find the truelocation from where a given particular read was derived, an alignmentfunction may be performed, e.g., a Smith-Waterman gapped alignment, aNeedleman-Wunsch alignment, etc., so as to determine where in the genomeone or more of the reads was actually derived, such as by comparing allof the possible locations where a match occurs and determining which ofall the possibilities is the most likely location in the genome fromwhich the read was sequenced, on the basis of which location's alignmentscore is greatest.

As indicated, typically, an algorithm is used to perform such analignment function. For example, a Smith-Waterman and/or aNeedleman-Wunsch alignment algorithm may be employed to align two ormore sequences against one another. In this instance, they may beemployed in a manner so as to determine the probabilities that for anygiven position where the read maps to the reference genome that themapping is in fact the position from where the read originated.Typically these algorithms are configured so as to be performed bysoftware, however, in various instances, such as herein presented, oneor more of these algorithms can be configured so as to be executed inhardware, as described in greater detail herein below.

In particular, the alignment function operates, at least in part, toalign one or more, e.g., all, of the reads to the reference genomedespite the presence of one or more portions of mismatches, e.g., SNPs,insertions, deletions, structural artifacts, etc. so as to determinewhere the reads are likely to fit in the genome correctly. For instance,the one or more reads are compared against the reference genome, and thebest possible fit for the read against the genome is determined, whileaccounting for substitutions and/or indels and/or structural variants.However, to better determine which of the modified versions of the readbest fits against the reference genome, the proposed changes must beaccounted for, and as such a scoring function may also be performed.

For instance, a scoring function may be performed, e.g., as part of anoverall alignment function, whereby as the alignment module performs itsfunction and introduces one or more changes into a sequence beingcompared to another, e.g., so as to achieve a better or best fit betweenthe two, for each change that is made so as to achieve the betteralignment, a number is detracted from a starting score, e.g., either aperfect score, or a zero starting score, in a manner such that as thealignment is performed the score for the alignment is also determined,such as where matches are detected the score is increased, and for eachchange introduced a penalty is incurred, and thus, the best fit for thepossible alignments can be determined, for example, by figuring outwhich of all the possible modified reads fits to the genome with thehighest score. Accordingly, in various instances, the alignment functionmay be configured to determine the best combination of changes that needto be made to the read(s) to achieve the highest scoring alignment,which alignment may then be determined to be the correct or most likelyalignment.

In view of the above, there are, therefore, at least two goals that maybe achieved from performing an alignment function. One is a report ofthe best alignment, including position in the reference genome and adescription of what changes are necessary to make the read match thereference segment at that position, and the other is the alignmentquality score. For instance, in various instances, the output from a thealignment module may be a Compact Idiosyncratic Gapped Alignment Report,e.g., a CIGAR string, wherein the CIGAR string output is a reportdetailing all the changes that were made to the reads so as to achievetheir best fit alignment, e.g., detailed alignment instructionsindicating how the query actually aligns with the reference. Such aCIGAR string readout may be useful in further stages of processing so asto better determine that for the given subject's genomic nucleotidesequence, the predicted variations as compared against a referencegenome are in fact true variations, and not just due to machine,software, or human error.

As set forth above, in various embodiments, alignment is typicallyperformed in a sequential manner, wherein the algorithm receives readsequence data, such as from a mapping module, pertaining to a read andone or more possible locations where the read may potentially map to theone or more reference genomes, and further receives genomic sequencedata, such as from one or more memories, pertaining to the one or morepositions in the one or more reference genomes to which the read maymap. In particular, in various embodiments, the mapping module processesthe reads, such as from a FASTQ file, and maps each of them to one ormore positions in the reference genome to where they may possibly align.The aligner then takes these predicted positions and uses them to alignthe reads to the reference genome, such as by building a virtual arrayby which the reads can be compared with the reference genome.

In performing this function the aligner evaluates each mapped positionfor each individual read and particularly evaluates those reads that mapto multiple possible locations in the reference genome and scores thepossibility that each position is the correct position. It then comparesthe best scores, e.g., the two best scores, and makes a decision as towhere the particular read actually aligns. For instance, in comparingthe first and second best alignment scores, the aligner looks at thedifference between the scores, and if the difference between them isgreat, then the confidence score that the one with the bigger score iscorrect will be high. However, where the difference between them issmall, e.g., zero, then the confidence score in being able to tell fromwhich of the two positions the read actually is derived is low, and moreprocessing may be useful in being able to clearly determine the truelocation in the reference genome from where the read is derived. Hence,the aligner in part is looking for the biggest difference between thefirst and second best confidence scores in making its call that a givenread maps to a given location in the reference genome. Ideally, thescore of the best possible choice of alignment is significantly greaterthan the score for the second best alignment for that sequence.

There are many different ways an alignment scoring methodology may beimplemented, for instance, each cell of the array may be scored or asub-portion of cells may be scored, such as in accordance with themethods disclosed herein. Typically, each alignment match, correspondingto a diagonal step in the alignment matrix, contributes a positivescore, such as +1, if the corresponding read and reference nucleotidesmatch; and a negative score, such as −4, if the two nucleotidesmismatch. Further, each deletion from the reference, corresponding to ahorizontal step in the alignment matrix, contributes a negative score,such as −7, and each insertion into the reference, corresponding to avertical step in the alignment matrix, contributes a negative score,such as −7.

In various instances, scoring parameters for nucleotide matches,nucleotide mismatches, insertions, and deletions may have any variouspositive or negative or zero values. In various instances, these scoringparameters may be modified based on available information. For instance,in certain instances, alignment gaps (insertions or deletions) arepenalized by an affine function of the gap length, for example −7 forthe first deleted (resp. inserted) nucleotide, but only −1 for eachadditional deleted (resp. inserted) nucleotide in continuous sequence.In various implementations, affine gap penalties may be achieved bysplitting gap (insertion or deletion) penalties into two components,such as a gap open penalty, e.g. −6, applied to the first step in a gap;and a gap extend penalty, e.g. −1, applied to every or further steps inthe gap. Affine gap penalties may yield more accurate alignments, suchas by letting alignments containing long insertions or deletions achieveappropriately high scores. Further, each lateral move may have the sameor different costs, such as the same cost per step, and/or where gapsoccur, such gaps can come at a higher or lower costs, such that the costfor lateral movements of the aligner may be less expensive than thecosts for gaps. Accordingly, in various embodiments, affine gap scoringmay be implemented, however, this can be expensive in software and/orhardware, because it typically requires a plurality, e.g., 3 scores, foreach cell to be scored, and hence, in various embodiments affine gapscoring is not implemented.

In various instances, scoring parameters may also be sensitive to “basequality scores” corresponding to nucleotides in the read. Some sequencedDNA read data, in formats such as FASTQ, may include a base qualityscore associated with each nucleotide, indicating an estimatedprobability that the nucleotide is incorrect, e.g. due to a sequencingerror. In some read data, base quality scores may indicate thelikelihood that an insertion and/or deletion sequencing error is presentin or adjacent to each position, or additional quality scores mayprovide this information separately. More accurate alignments,therefore, may be achieved by making scoring parameters, including anyor all of nucleotide match scores, nucleotide mismatch scores, gap(insertion and/or deletion) penalties, gap open penalties, and/or gapextend penalties, vary according to a base quality score associated withthe current read nucleotide or position. For example, score bonusesand/or penalties could be made smaller when a base quality scoreindicates a high probability a sequencing or other error being present.Base quality sensitive scoring may be implemented, for example, using afixed or configurable lookup-table, accessed using a base quality score,which returns corresponding scoring parameters.

In a hardware implementation in an integrated circuit, such as an FPGAor ASIC, a scoring wave front may be implemented as a linear array ofscoring cells, such as 16 cells, or 32 cells, or 64 cells, or 128 cellsor the like. Each of the scoring cells may be built of digital logicelements in a wired configuration to compute alignment scores. Hence,for each step of the wave front, for instance, each clock cycle, or someother fixed or variable unit of time, each of the scoring cells, or aportion of the cells, computes the score or scores required for a newcell in the virtual alignment matrix. Notionally, the various scoringcells are considered to be in various positions in the alignment matrix,corresponding to a scoring wave front as discussed herein, e.g., along astraight line extending from bottom-left to top-right in the matrix. Asis well understood in the field of digital logic design, the physicalscoring cells and their comprised digital logic need not be physicallyarranged in like manner on the integrated circuit.

Accordingly, as the wave front takes steps to sweep through the virtualalignment matrix, the notional positions of the scoring cellscorrespondingly update each cell, for example, notionally “moving” astep to the right, or for example, a step downward in the alignmentmatrix. All scoring cells make the same relative notional movement,keeping the diagonal wave front arrangement intact. Each time the wavefront moves to a new position, e.g., with a vertical downward step, or ahorizontal rightward step in the matrix, the scoring cells arrive in newnotional positions, and compute alignment scores for the virtualalignment matrix cells they have entered.

In such an implementation, neighboring scoring cells in the linear arrayare coupled to communicate query (read) nucleotides, referencenucleotides, and previously calculated alignment scores. The nucleotidesof the reference window may be fed sequentially into one end of the wavefront, e.g., the top-right scoring cell in the linear array, and mayshift from there sequentially down the length of the wave front, so thatat any given time, a segment of reference nucleotides equal in length tothe number of scoring cells is present within the cells, one successivenucleotide in each successive scoring cell.

Accordingly, each time the wave front steps horizontally, anotherreference nucleotide is fed into the top-right cell, and other referencenucleotides shift down-left through the wave front. This shifting ofreference nucleotides may be the underlying reality of the notionalmovement of the wave front of scoring cells rightward through thealignment matrix. Hence, the nucleotides of the read may be fedsequentially into the opposite end of the wave front, e.g. thebottom-left scoring cell in the linear array, and shift from theresequentially up the length of the wave front, so that at any given time,a segment of query nucleotides equal in length to the number of scoringcells is present within the cells, one successive nucleotide in eachsuccessive scoring cell.

Likewise, each time the wave front steps vertically, another querynucleotide is fed into the bottom-left cell, and other query nucleotidesshift up-right through the wave front.

This shifting of query nucleotides is the underlying reality of thenotional movement of the wave front of scoring cells downward throughthe alignment matrix. Accordingly, by commanding a shift of referencenucleotides, the wave front may be moved a step horizontally, and bycommanding a shift of query nucleotides, the wave front may be moved astep vertically. Accordingly, to produce generally diagonal wave frontmovement, such as to follow a typical alignment of query and referencesequences without insertions or deletions, wave front steps may becommanded in alternating vertical and horizontal directions.

Accordingly, neighboring scoring cells in the linear array may becoupled to communicate previously calculated alignment scores. Invarious alignment scoring algorithms, such as a Smith-Waterman orNeedleman-Wunsch, or such variant, the alignment score(s) in each cellof the virtual alignment matrix may be calculated using previouslycalculated scores in other cells of the matrix, such as the three cellspositioned immediately to the left of the current cell, above thecurrent cell, and diagonally up-left of the current cell. When a scoringcell calculates new score(s) for another matrix position it has entered,it must retrieve such previously calculated scores corresponding to suchother matrix positions. These previously calculated scores may beobtained from storage of previously calculated scores within the samecell, and/or from storage of previously calculated scores in the one ortwo neighboring scoring cells in the linear array. This is because thethree contributing score positions in the virtual alignment matrix(immediately left, above, and diagonally up-left) would have been scoredeither by the current scoring cell, or by one of its neighboring scoringcells in the linear array.

For instance, the cell immediately to the left in the matrix would havebeen scored by the current scoring cell, if the most recent wave frontstep was horizontal (rightward), or would have been scored by theneighboring cell down-left in the linear array, if the most recent wavefront step was vertical (downward). Similarly, the cell immediatelyabove in the matrix would have been scored by the current scoring cell,if the most recent wave front step was vertical (downward), or wouldhave been scored by the neighboring cell up-right in the linear array,if the most recent wave front step was horizontal (rightward).Similarly, the cell diagonally up-left in the matrix would have beenscored by the current scoring cell, if the most recent two wave frontsteps were in different directions, e.g., down then right, or right thendown, or would have been scored by the neighboring cell up-right in thelinear array, if the most recent two wave front steps were bothhorizontal (rightward), or would have been scored by the neighboringcell down-left in the linear array, if the most recent two wave frontsteps were both vertical (downward).

Accordingly, by considering information on the last one or two wavefront step directions, a scoring cell may select the appropriatepreviously calculated scores, accessing them within itself, and/orwithin neighboring scoring cells, utilizing the coupling betweenneighboring cells. In a variation, scoring cells at the two ends of thewave front may have their outward score inputs hard-wired to invalid, orzero, or minimum-value scores, so that they will not affect new scorecalculations in these extreme cells.

A wave front being thus implemented in a linear array of scoring cells,with such coupling for shifting reference and query nucleotides throughthe array in opposing directions, in order to notionally move the wavefront in vertical and horizontal steps, and coupling for accessingscores previously computed by neighboring cells in order to computealignment score(s) in new virtual matrix cell positions entered by thewave front, it is accordingly possible to score a band of cells in thevirtual matrix, the width of the wave front, such as by commandingsuccessive steps of the wave front to sweep it through the matrix. For anew read and reference window to be aligned, therefore, the wave frontmay begin positioned inside the scoring matrix, or, advantageously, maygradually enter the scoring matrix from outside, beginning e.g., to theleft, or above, or diagonally left and above the top-left corner of thematrix.

For instance, the wave front may begin with its top-left scoring cellpositioned just left of the top-left cell of the virtual matrix, and thewave front may then sweep rightward into the matrix by a series ofhorizontal steps, scoring a horizontal band of cells in the top-leftregion of the matrix. When the wave front reaches a predicted alignmentrelationship between the reference and query, or when matching isdetected from increasing alignment scores, the wave front may begin tosweep diagonally down-right, by alternating vertical and horizontalsteps, scoring a diagonal band of cells through the middle of thematrix. When the bottom-left wave front scoring cell reaches the bottomof the alignment matrix, the wave front may begin sweeping rightwardagain by successive horizontal steps, until some or all wave front cellssweep out of the boundaries of the alignment matrix, scoring ahorizontal band of cells in the bottom-right region of the matrix.

In a variation, increased efficiency may be obtained from the alignmentwave front by sharing its scoring cells between two successive alignmentoperations. A next alignment matrix having been established in advance,as the top-right portion of the wave front exits the bottom-right regionof the current alignment matrix, it may enter, immediately, or aftercrossing a minimum gap such as one cell or three cells, the top-rightregion of the next alignment matrix. In this manner, the horizontal wavefront sweep out of one alignment matrix can be the same motion as thehorizontal wave front sweep into the next alignment matrix. Doing thismay include the reference and query bases of the next alignment to befed into those scoring cells crossing into the next alignment matrix,and can reduce the average time consumed per alignment by the time toexecute a number of wave front steps almost equal to the number ofalignment cells in the wave front, e.g., such as 64 or 63 or 61 steps,which may take e.g. 64 or 63 or 61 clock cycles.

The number of scoring cells in an implementation of an alignment wavefront may be selected to balance various factors, including alignmentaccuracy, maximum insertion and deletion length, area, cost, and powerconsumption of the digital logic, clock frequency of the aligner logic,and performance of the overall integrated circuit. A long wave front isdesirable for good alignment accuracy, especially because a wave frontof N cells can align across indels approximately N nucleotides long, orslightly shorter. But a longer wave front costs more logic, whichconsumes more power. Further, a longer wave front can increase wirerouting complexity and delays on the integrated circuit, leading tolower maximum clock frequencies, reducing net aligner performance.Further still, if an integrated circuit has a limited size or powerconsumption, using a longer wave front may require less logic to beimplemented on the IC elsewhere, such as replicating fewer entire wavefronts, or other aligner or mapper logic components, this decreasing netperformance of the IC. In one particular embodiment, 64 scoring cells inthe wave front may give an acceptable balance of these factors.

Accordingly, where the wave front is X, e.g., 64 scoring cells wide, thescored band in the alignment matrix will likewise be 64 cells wide(measured diagonally). The matrix cells outside of this band do notnecessarily need to be processed nor their scores calculated, providedthat the optimal (best-scoring) alignment path through the matrix stayswithin the scored band. In a relatively small matrix, therefore, used toalign relatively short reads, e.g., 100 nucleotide or 250 nucleotidereads, this may be a safe assumption, such as if the wave front sweeps aperfect diagonal along the predicted aligned position of the read.

However, in some instances, such as in a large alignment matrix used toalign long reads, e.g., 1000 or 10,000 or 100,000 nucleotides, there maybe a substantial risk of accumulated indels causing the true alignmentto deviate from a perfect diagonal, sufficiently far in aggregate thatit may escape the scored band. In such instances, it may be useful tosteer the wave front so that the highest set of scores will be near thecenter of the wave front. Consequently, as the wave front performs itssweep, if the highest scores start to move one way or the other, e.g.,left to right, the wave front is shifted over to track this move. Forinstance, if the highest scores are observed in scoring cellssubstantially up-right from the center of the wave front, the wave frontmay be steered some distance straight rightward by successive horizontalsteps, until the highest scores return near the center of the wavefront.

Accordingly, an automatic steering mechanism may be implemented in thewave front control logic, to determine a steering target position withinthe length of the wave front, based on current and past scores observedin the wave front scoring cells, and to steer the wave front toward thistarget if it is off-center. More particularly, the position of themaximum score in the most recently scored wave front position may beused as a steering target. This is an effective method in someinstances. In some instances, however, the maximum score position may bea poor steering target. For instance, with some combinations ofalignment scoring parameters, when a long indel commences, and scoresaccordingly begin to decline, a pattern of two higher-score peaks with alower-score valley between them can form along the wave front, the twopeaks drifting apart as the indel continues.

Because it cannot be easily determined whether the event in progress isan insertion or a deletion, it is important for the wave front to trackdiagonally until successful matching commences again, either somedistance to the right for a deletion, or some distance downward for aninsertion. But if two spreading score peaks form, one of them is likelyto be slightly higher than the other, and could pull the automaticsteering in that direction, causing the wave front to lose the alignmentif the actual indel was in the other direction. A more robust method,therefore, may be to subtract a delta value from the maximum observedwave front score to determine a threshold score, identify the twoextreme scoring cells at least equal to this threshold score, and usethe midpoint between these extreme cells as the steering target. Thiswill tend to guide diagonally between a two-peak score pattern. Othersteering criteria can readily be applied, however, which serve to keephigher scores near the center of the wave front. If there is a delayedreaction between obtaining scores from wave front scoring cells andmaking a corresponding steering decision, hysteresis can advantageouslybe applied to compensate for steering decisions made in the interveningtime, to avoid oscillating patterns of automatic wave front steering.

One or more of such alignment procedures may be performed by anysuitable alignment algorithm, such as a Needleman-Wunsch alignmentalgorithm and/or a Smith-Waterman alignment algorithm that may have beenmodified to accommodate the functionality herein described. In generalboth of these algorithms and those like them basically perform, in someinstances, in a similar manner. For instance, as set forth above, thesealignment algorithms typically build the virtual array in a similarmanner such that, in various instances, the horizontal top boundary maybe configured to represent the genomic reference sequence, which may belaid out across the top row of the array according to its base paircomposition. Likewise, the vertical boundary may be configured torepresent the sequenced and mapped query sequences that have beenpositioned in order, downwards along the first column, such that theirnucleotide sequence order is generally matched to the nucleotidesequence of the reference to which they mapped. The intervening cellsmay then be populated with scores as to the probability that therelevant base of the query at a given position, is positioned at thatlocation relative to the reference. In performing this function, a swathmay be moved diagonally across the matrix populating scores within theintervening cells and the probability for each base of the query beingin the indicated position may be determined.

With respect to a Needleman-Wunsch alignment function, which generatesoptimal global (or semi-global) alignments, aligning the entire readsequence to some segment of the reference genome, the wave frontsteering may be configured such that it typically sweeps all the wayfrom the top edge of the alignment matrix to the bottom edge. When thewave front sweep is complete, the maximum score on the bottom edge ofthe alignment matrix (corresponding to the end of the read) is selected,and the alignment is back-traced to a cell on the top edge of the matrix(corresponding to the beginning of the read). In various of theinstances disclosed herein, the reads can be any length long, can be anysize, and there need not be extensive read parameters as to how thealignment is performed, e.g., in various instances, the read can be aslong as a chromosome. In such an instance, however, the memory size andchromosome length may be limiting factor.

With respect to a Smith-Waterman algorithm, which generates optimallocal alignments, aligning the entire read sequence or part of the readsequence to some segment of the reference genome, this algorithm may beconfigured for finding the best scoring possible based on a full orpartial alignment of the read. Hence, in various instances, the wavefront-scored band may not extend to the top and/or bottom edges of thealignment matrix, such as if a very long read had only seeds in itsmiddle mapping to the reference genome, but commonly the wave front maystill score from top to bottom of the matrix. Local alignment istypically achieved by two adjustments. First, alignment scores are neverallowed to fall below zero (or some other floor), and if a cell scoreotherwise calculated would be negative, a zero score is substituted,representing the start of a new alignment. Second, the maximum alignmentscore produced in any cell in the matrix, not necessarily along thebottom edge, is used as the terminus of the alignment. The alignment isbacktraced from this maximum score up and left through the matrix to azero score, which is used as the start position of the local alignment,even if it is not on the top row of the matrix.

In view of the above, there are several different possible pathwaysthrough the virtual array. In various embodiments, the wave front startsfrom the upper left corner of the virtual array, and moves downwardstowards identifiers of the maximum score. For instance, the results ofall possible aligns can be gathered, processed, correlated, and scoredto determine the maximum score. When the end of a boundary or the end ofthe array has been reached and/or a computation leading to the highestscore for all of the processed cells is determined (e.g., the overallhighest score identified) then a backtrace may be performed so as tofind the pathway that was taken to achieve that highest score.

For example, a pathway that leads to a predicted maximum score may beidentified, and once identified an audit may be performed so as todetermine how that maximum score was derived, for instance, by movingbackwards following the best score alignment arrows retracing thepathway that led to achieving the identified maximum score, such ascalculated by the wave front scoring cells. This backwardsreconstruction or backtrace involves starting from a determined maximumscore, and working backward through the previous cells navigating thepath of cells having the scores that led to achieving the maximum scoreall the way up the table and back to an initial boundary, such as thebeginning of the array, or a zero score in the case of local alignment.

During a backtrace, having reached a particular cell in the alignmentmatrix, the next backtrace step is to the neighboring cell, immediatelyleftward, or above, or diagonally up-left, which contributed the bestscore that was selected to construct the score in the current cell. Inthis manner, the evolution of the maximum score may be determined,thereby figuring out how the maximum score was achieved. The backtracemay end at a corner, or an edge, or a boundary, or may end at a zeroscore, such as in the upper left hand corner of the array. Accordingly,it is such a back trace that identifies the proper alignment and therebyproduces the CIGAR strand readout, e.g., 3M, 2D, 8M, 4I, 16M, etc., thatrepresents how the sample genomic sequence derived from the individual,or a portion thereof, matches to, or otherwise aligns with, the genomicsequence of the reference DNA.

Accordingly, once it has been determined where each read is mapped, andfurther determined where each read is aligned, e.g., each relevant readhas been given a position and a quality score reflecting the probabilitythat the position is the correct alignment, such that the nucleotidesequence for the subject's DNA is known, then the order of the variousreads and/or genomic nucleic acid sequence of the subject may beverified, such as by performing a back trace function moving backwardsup through the array so as to determine the identity of every nucleicacid in its proper order in the sample genomic sequence. Consequently,in some aspects, the present disclosure is directed to a back tracefunction, such as is part of an alignment module that performs both analignment and a back trace function, such as a module that may be partof a pipeline of modules, such as a pipeline that is directed at takingraw sequence read data, such as form a genomic sample form anindividual, and mapping and/or aligning that data, which data may thenbe sorted.

To facilitate the backtrace operation, it is useful to store a scoringvector for each scored cell in the alignment matrix, encoding thescore-selection decision. For classical Smith-Waterman and/orNeedleman-Wunsch scoring with linear gap penalties, the scoring vectorcan encode four possibilities, which may optionally be stored as a 2-bitinteger from 0 to 3, for example: 0=new alignment (null score selected);1=vertical alignment (score from the cell above selected, modified bygap penalty); 2=horizontal alignment (score from the cell to the leftselected, modified by gap penalty); 3=diagonal alignment (score from thecell up and left selected, modified by nucleotide match or mismatchscore). Optionally, the computed score(s) for each scored matrix cellmay also be stored (in addition to the maximum achieved alignment scorewhich is standardly stored), but this is not generally necessary forbacktrace, and can consume large amounts of memory. Performing backtracethen becomes a matter of following the scoring vectors; when thebacktrace has reached a given cell in the matrix, the next backtracestep is determined by the stored scoring vector for that cell, e.g.:0=terminate backtrace; 1=backtrace upward; 2=backtrace leftward;3=backtrace diagonally up-left.

Such scoring vectors may be stored in a two-dimensional table arrangedaccording to the dimensions of the alignment matrix, wherein onlyentries corresponding to cells scored by the wave front are populated.Alternatively, to conserve memory, more easily record scoring vectors asthey are generated, and more easily accommodate alignment matrices ofvarious sizes, scoring vectors may be stored in a table with each rowsized to store scoring vectors from a single wave front of scoringcells, e.g. 128 bits to store 64 2-bit scoring vectors from a 64-cellwave front, and a number of rows equal to the maximum number of wavefront steps in an alignment operation.

Additionally, for this option, a record may be kept of the directions ofthe various wavefront steps, e.g., storing an extra, e.g., 129^(th), bitin each table row, encoding e.g., 0 for vertical wavefront steppreceding this wavefront position, and 1 for horizontal wavefront steppreceding this wavefront position. This extra bit can be used duringbacktrace to keep track of which virtual scoring matrix positions thescoring vectors in each table row correspond to, so that the properscoring vector can be retrieved after each successive backtrace step.When a backtrace step is vertical or horizontal, the next scoring vectorshould be retrieved from the previous table row, but when a backtracestep is diagonal, the next scoring vector should be retrieved from tworows previous, because the wavefront had to take two steps to move fromscoring any one cell to scoring the cell diagonally right-down from it.

In the case of affine gap scoring, scoring vector information may beextended, e.g. to 4 bits per scored cell. In addition to the e.g. 2-bitscore-choice direction indicator, two 1-bit flags may be added, avertical extend flag, and a horizontal extend flag. According to themethods of affine gap scoring extensions to Smith-Waterman orNeedleman-Wunsch or similar alignment algorithms, for each cell, inaddition to the primary alignment score representing the best-scoringalignment terminating in that cell, a ‘vertical score’ should begenerated, corresponding to the maximum alignment score reaching thatcell with a final vertical step, and a ‘horizontal score’ should begenerated, corresponding to the maximum alignment score reaching thatcell with a final horizontal step; and when computing any of the threescores, a vertical step into the cell may be computed either using theprimary score from the cell above minus a gap-open penalty, or using thevertical score from the cell above minus a gap-extend penalty, whicheveris greater; and a horizontal step into the cell may be computed eitherusing the primary score from the cell to the left minus a gap-openpenalty, or using the horizontal score from the cell to the left minus agap-extend penalty, whichever is greater. In cases where the verticalscore minus a gap extend penalty is selected, the vertical extend flagin the scoring vector should be set, e.g. ‘1’, and otherwise it shouldbe unset, e.g. ‘0’. In cases when the horizontal score minus a gapextend penalty is selected, the horizontal extend flag in the scoringvector should be set, e.g. ‘1’, and otherwise it should be unset, e.g.‘0’. During backtrace for affine gap scoring, any time backtrace takes avertical step upward from a given cell, if that cell's scoring vector'svertical extend flag is set, the following backtrace step must also bevertical, regardless of the scoring vector for the cell above. Likewise,any time backtrace takes a horizontal step leftward from a given cell,if that cell's scoring vector's horizontal extend flag is set, thefollowing backtrace step must also be horizontal, regardless of thescoring vector for the cell to the left.

Accordingly, such a table of scoring vectors, e.g. 129 bits per row for64 cells using linear gap scoring, or 257 bits per row for 64 cellsusing affine gap scoring, with some number NR of rows, is adequate tosupport backtrace after concluding alignment scoring where the scoringwavefront took NR steps or fewer. For example, when aligning300-nucleotide reads, the number of wavefront steps required may alwaysbe less than 1024, so the table may be 257×1024 bits, or approximately32 kilobytes, which in many cases may be a reasonable local memoryinside the IC. But if very long reads are to be aligned, e.g. 100,000nucleotides, the memory requirements for scoring vectors may be quitelarge, e.g. 8 megabytes, which may be very costly to include as localmemory inside the IC. For such support, scoring vector information maybe recorded to bulk memory outside the IC, e.g. DRAM, but then thebandwidth requirements, e.g. 257 bits per clock cycle per alignermodule, may be excessive, which may bottleneck and dramatically reducealigner performance.

Accordingly, it is desirable to have a method for disposing of scoringvectors before completing alignment, so their storage requirements canbe kept bounded, e.g. to perform incremental backtraces, generatingincremental partial CIGAR strings for example, from early portions of analignment's scoring vector history, so that such early portions of thescoring vectors may then be discarded. The challenge is that thebacktrace is supposed to begin in the alignment's terminal, maximumscoring cell, which unknown until the alignment scoring completes, soany backtrace begun before alignment completes may begin from the wrongcell, not along the eventual final optimal alignment path.

Accordingly, a method is given for performing incremental backtrace frompartial alignment information, e.g. comprising partial scoring vectorinformation for alignment matrix cells scored so far. From a currentlycompleted alignment boundary, e.g., a particular scored wave frontposition, backtrace is initiated from all cell positions on theboundary. Such backtrace from all boundary cells may be performedsequentially, or advantageously, especially in a hardwareimplementation, all the backtraces may be performed together. It is notnecessary to extract alignment notations, e.g., CIGAR strings, fromthese multiple backtraces; only to determine what alignment matrixpositions they pass through during the backtrace. In an implementationof simultaneous backtrace from a scoring boundary, a number of 1-bitregisters may be utilized, corresponding to the number of alignmentcells, initialized e.g., all to ‘1’s, representing whether any of thebacktraces pass through a corresponding position. For each step ofsimultaneous backtrace, scoring vectors corresponding to all the current‘1’s in these registers, e.g. from one row of the scoring vector table,can be examined, to determine a next backtrace step corresponding toeach ‘1’ in the registers, leading to a following position for each ‘1’in the registers, for the next simultaneous backtrace step.

Importantly, it is easily possible for multiple 1's in the registers tomerge into common positions, corresponding to multiple of thesimultaneous backtraces merging together onto common backtrace paths.Once two or more of the simultaneous backtraces merge together, theyremain merged indefinitely, because henceforth they will utilize scoringvector information from the same cell. It has been observed, empiricallyand for theoretical reasons, that with high probability, all of thesimultaneous backtraces merge into a singular backtrace path, in arelatively small number of backtrace steps, which e.g. may be a smallmultiple, e.g. 8, times the number of scoring cells in the wavefront.For example, with a 64-cell wavefront, with high probability, allbacktraces from a given wavefront boundary merge into a single backtracepath within 512 backtrace steps. Alternatively, it is also possible, andnot uncommon, for all backtraces to terminate within the number, e.g.512, of backtrace steps.

Accordingly, the multiple simultaneous backtraces may be performed froma scoring boundary, e.g. a scored wavefront position, far enough backthat they all either terminate or merge into a single backtrace path,e.g. in 512 backtrace steps or fewer. If they all merge together into asingular backtrace path, then from the location in the scoring matrixwhere they merge, or any distance further back along the singularbacktrace path, an incremental backtrace from partial alignmentinformation is possible. Further backtrace from the merge point, or anydistance further back, is commenced, by normal singular backtracemethods, including recording the corresponding alignment notation, e.g.,a partial CIGAR string. This incremental backtrace, and e.g. partialCIGAR string, must be part of any possible final backtrace, and e.g.full CIGAR string, that would result after alignment completes, unlesssuch final backtrace would terminate before reaching the scoringboundary where simultaneous backtrace began, because if it reaches thescoring boundary, it must follow one of the simultaneous backtracepaths, and merge into the singular backtrace path, now incrementallyextracted.

Therefore, all scoring vectors for the matrix regions corresponding tothe incrementally extracted backtrace, e.g., in all table rows for wavefront positions preceding the start of the extracted singular backtrace,may be safely discarded. When the final backtrace is performed from amaximum scoring cell, if it terminates before reaching the scoringboundary (or alternatively, if it terminates before reaching the startof the extracted singular backtrace), the incremental alignmentnotation, e.g. partial CIGAR string, may be discarded. If the finalbacktrace continues to the start of the extracted singular backtrace,its alignment notation, e.g., CIGAR string, may then be grafted onto theincremental alignment notation, e.g., partial CIGAR string.

Furthermore, in a very long alignment, the process of performing asimultaneous backtrace from a scoring boundary, e.g., scored wave frontposition, until all backtraces terminate or merge, followed by asingular backtrace with alignment notation extraction, may be repeatedmultiple times, from various successive scoring boundaries. Theincremental alignment notation, e.g. partial CIGAR string, from eachsuccessive incremental backtrace may then be grafted onto theaccumulated previous alignment notations, unless the new simultaneousbacktrace or singular backtrace terminates early, in which caseaccumulated previous alignment notations may be discarded. The eventualfinal backtrace likewise grafts its alignment notation onto the mostrecent accumulated alignment notations, for a complete backtracedescription, e.g. CIGAR string.

Accordingly, in this manner, the memory to store scoring vectors may bekept bounded, assuming simultaneous backtraces always merge together ina bounded number of steps, e.g. 512 steps. In rare cases wheresimultaneous backtraces fail to merge or terminate in the bounded numberof steps, various exceptional actions may be taken, including failingthe current alignment, or repeating it with a higher bound or with nobound, perhaps by a different or traditional method, such as storing allscoring vectors for the complete alignment, such as in external DRAM. Ina variation, it may be reasonable to fail such an alignment, because itis extremely rare, and even rarer that such a failed alignment wouldhave been a best-scoring alignment to be used in alignment reporting.

In an optional variation, scoring vector storage may be divided,physically or logically, into a number of distinct blocks, e.g. 512 rowseach, and the final row in each block may be used as a scoring boundaryto commence a simultaneous backtrace. Optionally, a simultaneousbacktrace may be required to terminate or merge within the single block,e.g. 512 steps. Optionally, if simultaneous backtraces merge in fewersteps, the merged backtrace may nevertheless be continued through thewhole block, before commencing an extraction of a singular backtrace inthe previous block. Accordingly, after scoring vectors are fully writtento block N, and begin writing to block N+1, a simultaneous backtrace maycommence in block N, followed by a singular backtrace and alignmentnotation extraction in block N−1. If the speed of the simultaneousbacktrace, the singular backtrace, and alignment scoring are all similaror identical, and can be performed simultaneously, e.g., in parallelhardware in an IC, then the singular backtrace in block N−1 may besimultaneous with scoring vectors filling block N+2, and when block N+3is to be filled, block N−1 may be released and recycled.

Thus, in such an implementation, a minimum of 4 scoring vector blocksmay be employed, and may be utilized cyclically. Hence, the totalscoring vector storage for an aligner module may be 4 blocks of 257×512bits each, for example, or approximately 64 kilobytes. In a variation,if the current maximum alignment score corresponds to an earlier blockthan the current wavefront position, this block and the previous blockmay be preserved rather than recycled, so that a final backtrace maycommence from this position if it remains the maximum score; having anextra 2 blocks to keep preserved in this manner brings the minimum,e.g., to 6 blocks. In another variation, to support overlappedalignments, the scoring wave front crossing gradually from one alignmentmatrix to the next as described above, additional blocks, e.g. 1 or 2additional blocks, may be utilized, e.g., 8 blocks total, e.g.,approximately 128 kilobytes. Accordingly, if such a limited number ofblocks, e.g., 4 blocks or 8 blocks, is used cyclically, alignment andbacktrace of arbitrarily long reads is possible, e.g., 100,000nucleotides, or an entire chromosome, without the use of external memoryfor scoring vectors.

It is to be understood, such as with reference to the above, thatalthough a mapping function may in some instances have been described,such as with reference to a mapper, and/or an alignment function mayhave in some instances been described, such as with reference to analigner, these different functions may be performed sequentially by thesame architecture, which has commonly been referenced in the art as analigner. Accordingly, in various instances, both the mapping functionand the aligning function, as herein described may be performed by acommon architecture that may be understood to be an aligner, especiallyin those instances wherein to perform an alignment function, a mappingfunction need first be performed.

The output from the alignment module is a SAM (Text) or BAM (e.g.,binary version of a SAM) file along with a mapping quality score (MAPA),which quality score reflects the confidence that the predicted andaligned location of the read to the reference is actually where the readis derived. Accordingly, once it has been determined where each read ismapped, and further determined where each read is aligned, e.g., eachrelevant read has been given a position and a quality score reflectingthe probability that the position is the correct alignment, such thatthe nucleotide sequence for the subject's DNA is known as well as howthe subject's DNA differs from that of the reference (e.g., the CIGARstring has been determined), then the various reads representing thegenomic nucleic acid sequence of the subject may be sorted by chromosomelocation, so that the exact location of the read on the chromosomes maybe determined. Consequently, in some aspects, the present disclosure isdirected to a sorting function, such as may be performed by a sortingmodule, which sorting module may be part of a pipeline of modules, suchas a pipeline that is directed at taking raw sequence read data, such asform a genomic sample form an individual, and mapping and/or aligningthat data, which data may then be sorted.

More particularly, once the reads have been assigned a position, such asrelative to the reference genome, which may include identifying to whichchromosome the read belongs and/or its offset from the beginning of thatchromosome, the reads may be sorted by position. Sorting may be useful,such as in downstream analyses, whereby all of the reads that overlap agiven position in the genome may be formed into a pile up so as to beadjacent to one another, such as after being processed through thesorting module, whereby it can be readily determined if the majority ofthe reads agree with the reference value or not. Hence, where themajority of reads do not agree with the reference value a variant callcan be flagged. Sorting, therefore, may involve one or more of sortingthe reads that align to the relatively same position, such as the samechromosome position, so as to produce a pileup, such that all the readsthat cover the same location are physically grouped together; and mayfurther involve analyzing the reads of the pileup to determine where thereads may indicate an actual variant in the genome, as compared to thereference genome, which variant may be distinguishable, such as by theconsensus of the pileup, from an error, such as a machine read error orerror an error in the sequencing methods which may be exhibited by asmall minority of the reads.

Once the data has been obtained there are one or more other modules thatmay be run so as to clean up the data. For instance, one module that maybe included, for example, in a sequence analysis pipeline, such as fordetermining the genomic sequence of an individual, may be a localrealignment module. For example, it is often difficult to determineinsertions and deletions that occur at the end of the read. This isbecause the Smith-Waterman or equivalent alignment process lacks enoughcontext beyond the indel to allow the scoring to detect its presence.Consequently, the actual indel may be reported as one or more SNPs. Insuch an instance, the accuracy of the predicted location for any givenread may be enhanced by performing a local realignment on the mappedand/or aligned and/or sorted read data.

In such instances, pileups may be used to help clarify the properalignment, such as where a position in question is at the end of anygiven read, that same position is likely to be at the middle of someother read in the pileup. Accordingly, in performing a local realignmentthe various reads in a pileup may be analyzed so as to determine if someof the reads in the pile up indicate that there was an insertion or adeletion at a given position where an other read does not include theindel, or rather includes a substitution, at that position, then theindel may be inserted, such as into the reference, where it is notpresent, and the reads in the local pileup that overlap that region maybe realigned to see if collectively a better score is achieved then whenthe insertion and/or deletion was not there. Accordingly, if there is animprovement, the whole set of reads in the pileup may be reviewed and ifthe score of the overall set has improved then it is clear to make thecall that there really was an indel at that position. In a manner suchas this, the fact that there is not enough context to more accuratelyalign a read at the end of a chromosome, for any individual read, may becompensated for. Hence, when performing a local realignment, one or morepileups where one or more indels may be positioned are examined, and itis determined if by adding an indel at any given position the overallalignment score may be enhanced.

Another module that may be included, for example, in a sequence analysispipeline, such as for determining the genomic sequence of an individual,may be a duplicate marking module. For instance, a duplicate markingfunction may be performed so as to compensate for chemistry errors thatmay occur during the sequencing phase. For example, as described above,during some sequencing procedures nucleic acid sequences are attached tobeads and built up from there using labeled nucleotide bases. Ideallythere will be only one read per bead. However, sometimes multiple readsbecome attached to a single bead and this results in an excessive numberof copies of the attached read. This phenomenon is known as readduplication.

Such read duplication may throw off the statistics and create astatistical bias because instead of having an equal representation ofall reads, various reads have been duplicated, such as because of theduplicate template sequences attached to more than one bead are overrepresented. Accordingly, these may be determined because any read thataligns to the exact same position, and has the exact same length, islikely a duplicate. Once this is identified by the system, only one readneed be subjected to further processing and the others may be marked asduplicates and, therefore, can be discarded or ignored. A typicalsituation where this occurs is where there is not enough geneticmaterial to process from the very beginning and the system attempts toovercompensate for that.

Another module that may be included, for example, in a sequence analysispipeline, such as for determining the genomic sequence of an individual,may be a base quality score recalibrater. For instance, every base ofevery read has a Phred score that indicates the probability that thecalled base at that position is incorrect. For example, the Phred scorefor any base is due in part to the nature of the base that precedes itand the error profile will be different depending on which base precedesthe base in question. Further, there is a greater likelihood of an erroroccurring at the ends of a read, e.g., such as where at the ends of thereads the chemistry is starting to lose its performance. A base qualityscore recalibration is a covariant analysis that may go back andmeasures the empirical quality of the base quality score as a functionof all those things by which it varies.

In various instances, it involves two passes, the first gathers all theactual, empirical measured data and statistics on the error rateobserved as a function of all the variables, and the second passinvolves the actual recalibration of the scores by flowing all the readsthrough a filter modifying the quality scores for every single base as afunction of the variables based on what was actually empiricallymeasured in the data set. This compensates for all the differences inthe data due to the various variables and cleans up that data and score.The purpose of all this cleanup is to ensure the best possible variantcalling is achieved. Many variant callers base their decisions in parton the reported quality of each of the nucleotides that pile up at eachposition in the genome. If the quality scores are not accurate, therecould easily result a wrong call.

Another module that may be included, for example, in a sequence analysispipeline, such as for determining the genomic sequence of an individual,may be a compression module, that executes a compression function. Asindicated above, it may be useful at some point to take the generatedand processed data and transmit it to a remote location, such as thecloud, and hence, the data may need to be compressed at a particularstage of processing, whereby once compressed it may be transmittedand/or otherwise uploaded, such as on to the cloud or to a server farm,etc., for instance, for the performance of the variant calling module.The results once obtained may then be decompressed and/or stored in thememory, on a data base on the cloud, such as an electronic health and/orresearch database, and the like, which in turn, can be made availablefor tertiary processing, etc.

Accordingly, as set forth herein above, in various aspects, this presentdisclosure is directed to systems, apparatuses, and methods forimplementing genomics and/or bioinformatic protocols such as, in variousinstances, for performing one or more functions for analyzing geneticdata on an integrated circuit, such as implemented in a hardwareprocessing platform. For example, in one aspect, a bioinformatics systemis provided, wherein the system may involve the performance of variousbioanalytical functions that have been optimized so as to be performedfaster and/or with increased accuracy in a hardware implementation.Accordingly, in various instances, the methods and systems hereindescribed may include the performance of one or more algorithms forexecuting these functions, wherein the algorithms may be implemented ina hardware solution, such as where the algorithm has been optimized soas to be implemented by an integrated circuit formed of one or morehardwired digital logic circuits. In such an instance, the hardwireddigital logic circuits may be interconnected, such as by one or aplurality of physical electrical interconnects, and may be arranged tofunction as one or more processing engines. In various instances, aplurality of hardwired digital logic circuits are provided, whichhardwired digital logic circuits are configured as a set of processingengines, wherein each processing engine is capable of performing one ormore steps in the bioinformatics genetic analysis protocol.

More particularly, in one instance, a system for executing a sequenceanalysis pipeline such as on genetic sequence data is provided. Thesystem may include one or more of an electronic data source, a memory,and an integrated circuit. For instance, in one embodiment, anelectronic data source is included, where in the electronic data sourcemay be configured for providing one or more digital signals, such as adigital signal representing one or more reads of genetic data, forexample, where each read of genomic data includes a sequence ofnucleotides. Further, the memory may be configured for storing one ormore genetic reference sequences, and may further be configured forstoring an index, such as an index of the one or more genetic referencesequences.

Further still, in various instances, one or more of the plurality ofphysical electrical interconnects may include an input, such as to theintegrated circuit, and may further be connected with the electronicdata source, so as to be able to receive the one or more reads ofgenomic data. In various embodiments, the hardwired digital logiccircuits may be arranged as a set of processing engines, such as whereeach processing engine is formed of a subset of the hardwired digitallogic circuits, and is configured so as to perform one or more steps inthe sequence analysis pipeline, such as on digitized genetic data, e.g.,on the plurality of reads of genomic data. In such instances, eachsubset of the hardwired digital logic circuits may be in a wiredconfiguration so as to perform the one or more steps in the sequenceanalysis pipeline, such as where the one or more steps may includeperforming one or more of: a base calling and/or error correctionoperation, such as on the digitized genetic data, and/or may include oneor more of performing a mapping, an alignment, and/or a sorting functionon the genetic data. In certain instances, the pipeline may includeperforming one or more of a realignment, a deduplication, a base qualityscore recalibration, a reduction and/or compression, and/or adecompression on the digitized genetic data. In certain instances thepipeline may include performing a variant calling operation on thegenetic data.

Accordingly, in various embodiments, the systems, apparatuses, andmethods for implementing genomics and/or bioinformatic protocols, asherein described, may involve taking processes that may have typicallybeen performed on software, and embedding those functions into anintegrated circuit, such as on a chip, for instance as part of a circuitboard, such as where the functions have been optimized to enhance itsperformance on the chip. Hence, in one embodiment, as can be seen withrespect to FIG. 1 a chip is provided wherein the chip has been designedso as to efficiently perform the functions of the pipeline. In variousparticular embodiments the chip may be a field programmable gate array(FPGA), or an application specific integrated circuit (ASIC), or thelike.

For instance, the functioning of one or more of these algorithms may beembedded onto a chip, such as into an FPGA or ASIC chip, and may beoptimized so as to perform more efficiently because of theirimplementation in such hardware. Accordingly, in one embodiment a FPGAchip is provided wherein the chip is capable of being configurable,e.g., its programming may be changed, so as to be more adaptable inmeeting a given user's needs with respect to performing the variousgenomic functions detailed herein. In such an instance, the user canchange and/or modify the algorithms employed dependent on the keyparameters desired to be emphasized in the overall system, such as togive additional functionality or change out what was first presented onthe chip, e.g., such as re-configuring the chip to employ a differentalgorithm. In accordance with another embodiment an ASIC is provided,such as where the FPGA is converted to an ASIC chip where itsfunctionality is locked down into the chip. In such an instance, variousparameters, such as various parameters regarding the function of one ormore of the algorithms set forth herein, may be user selected, forinstance, governing how the various modules are supposed to function,but the way those modules actually function is locked in.

In various embodiments, as seen with respect to FIG. 1, the chip may bepart of a circuit board, such as part of an expansion card, forinstance, a peripheral component interconnect (PCI) card, including aPCIe card, which in various embodiments may be associated, such as,communicably coupled, e.g., electrically connected, with an automatedsequencer device so as to function part and parcel with the sequencer,such as where the data files, e.g., FASTQ files, generated by thesequencer is transferred directly over to the chip, such as forsecondary genomic processing, such as immediately subsequent to theFASTQ file generation and/or primary processing, e.g., immediately afterthe sequencing function has been performed.

Accordingly, in certain instances, a PCI card is provided wherein thePCI card may include a chip with a PCIe bus, where the chip may includeone or more of: a configuration manager, such as a configuration control(Cent-Com); a direct memory access engine (e.g., a driver); an API; aclient level interface (CLI), a library; a memory, such as a randomaccess memory (RAM) or a dynamic random access memory (DRAM); and/or achip level interconnect, such as a DDR3. For instance, in variousinstances a configuration manager may be included wherein theconfiguration manager is driven, such as by a parameter file. In such aninstance the configuration manager may be adapted so as to configure thevarious modules of the pipeline. In various instances, it may be usereditable, and thereby allow a user to determine which modules of thepipeline are going to be used, e.g., from all of them to a subset ofless than all of them, such as for a particular dataset, such as aparticular set of FASTQ files.

For example, in various embodiments, the functioning of the pipeline isvery configurable such that one or more of the modules, such asstructured into the chip, may be run or not run, as desired. Further,each module in use can also be configured so as to run in accordancewith one or more preselected parameters, which the user may have controlover, such as regarding how the module is going to perform and behave.Hence, there may be two different sets of configuration files, such asone that controls the basic operations of the system as a whole, and maybe hidden from the user, and another that is capable of beingmanipulated by the user, thereby allowing the user to select various ofthe parameters by which one or more of the subsystems, e.g., modules, ofthe chip will be run.

Further still, various of the above described modules may be hardwiredinto the chip, or may be external to the chip, but positioned in acoupling relationship therewith, such as on a PCI board, or they may belocated remotely from the chip, such as on a different PCI board, oreven on a different server, such as on a server that may be accessed viathe cloud. For instance, in certain implementations, one or more of theabove described modules may be hardwired onto a chip and the chipinstalled onto the circuit board of a stand-alone device, or coupled toa sequencer, whereby the user configures and runs the system directly bythemselves according to their own preselected parameters. Alternatively,as indicated herein, one or more of the above described modules may bepresent on a system that is accessible via the cloud, wherein thedirecting of the functioning of the pipeline, and/or the modulesthereof, may include the user logging on to a server, e.g., a remoteserver, and transmitting data to and therefrom, and thereby selectswhich modules to be run on the data set. In certain instances, one ormore of the modules may be performed remotely, such as via the cloudaccessed server.

In various instances, in configuring the system, the chip, e.g., thechip on an expansion card, such as a PCI card, may be included in aserver, whereby the server runs the various applications of the system.In certain instances, the server may have a terminal connectable therewith, whereby a windows interface may be presentable to the user suchthat the user may select the modules to be run and the parameters bywhich they are to be run, such as by selecting a box from a menu ofboxes. In other instances, however, the parameter file may be a textfile detailing categories by module under file names that the user canthen edit, so as to select which modules will be run in accordance withwhich parameters. For instance, in various embodiments, each chip mayinclude all or a selection of the modules, such as one or more of: abase calling, error correcting, a mapping, an alignment, a sorting, alocal realignment, a duplicate marking, a recalibration, a variantcalling, a compression, and/or a decompression module, from which theuser may select which modules will run, when, and to various extents howit will run, without changing the functioning of the underlyingalgorithms by which the individual modules are operated.

Additionally, in various instances, a direct memory access (DMA) enginein the chip, and a DMA driver, may be included wherein the DMA driverincludes code that runs in the kernel. Accordingly, the DMA driver maybe the foundation of the overall operating system. For instance, wherethe kernel runs in a literal addressing space, layered above that may bea virtual user space. This operating system software, therefore operatesin between these layers managing the mapping from the virtual to thephysical space. More particularly, the kernel represents the lowestlevel of code that gives the platform access to the PCI, e.g., PCIe,bus, to which the chip is coupled. Accordingly, since, in variousembodiments, the chip may be configured as an expansion card with a PCIeexpansion bus, which expansion card may be coupled with various hardwareof a device, such as a sequencer, the DMA driver may function so as tocommunicate with the hardware of the sequencer, and may further beconfigured for running at the kernel level on the CPU, so as to alsocommunicate with the DMA engine in the chip, and/or be configured foroperating in the virtual user space so as to receive instructions fromthe user.

To facilitate this communication within the chip and/or between the chipand one or more cards, every single configurable parameter of a modulemay be assigned to a register address. In such an instance, the card mayhave its own address space, which address space may be different fromthe address space for one or more memories, such as 64 gigabytes ofmemory, and/or additionally every module may have registers and localmemory associated with it, each with its own address space. Accordingly,the driver knows where everything is, all the addresses, and knows howto communicate between the chip, the PCI card, and/or the hardware ofthe server. Further, knowing where all the addresses are andcommunicating with an API the driver can read the parameter file that auser generates, and can look up for that parameter where the file isactually located in the host computer system and will read and interpretthe value in the file and will deliver that value in the right registerin the right place in the chip. Hence, the driver may handle deliveringthe selected parameter instructions, such as with respect to varioususer selected configurations, and ships that data to the chip via theDMA engine to configure any of its processing functions.

Further, in various instances, an API may be included wherein the API isconfigured so as to include a list of function calls that the user canmake, so as to configure and operate the system. For instance, an APImay be defined in a header file that describes the functionality anddetermines how to call a function, such as the parameters that arepassed, the inputs and outputs, what comes in, what goes out, and whatgets returned. For example, in various embodiments, one or more of theelements of the pipeline may be configurable such as by instructionsentered by a user and/or one or more third party applications. Theseinstructions may be communicated to the chip via the API whichcommunicates with the driver, instructing the driver as to which partsof the chip, e.g., which modules are to be activated, when, and in whatorder, given a preselected parameter configuration.

As indicated above, the DMA driver runs at the kernel level, and has itsown very low level, basic API that provides access to the hardware andfunctions so as to access applicable registers and modules. On top ofthis layer is built a virtual layer of service functions, that form thebuilding blocks that are used for a multiplicity of functions that sendfiles down to the kernel and gets results back, and further performsmore higher level functions. On top of that layer is an additional layerthat uses those service functions, which is the API level that a userwill interface with and it functions primarily for configuration,downloading files, and uploading results. Such configuration may includecommunicating with registers and also performing function calls.

For example, as described herein above, one function call may be togenerate the hash table via the hashing algorithm. Specifically, becausein certain embodiments this function may be based on a reference genome,once for every reference genome, the hash tables that are used in themapper may need to be constructed, based on the reference, there istherefore a function call that performs this function, which functioncall will accept a file name of where the reference file is stored andit will then generate one or more data files that contain the hash tableand the reference. Another function call may be to load the hash tablethat was generated via the hashing algorithm and transfer that down tothe memory on the chip, and/or put it at the right spot where thehardware is expecting them to be. Of course, the reference itself willneed to be downloaded onto the chip, as well for the performance of thealignment function, and the configuration manager can perform thatfunction such as by loading everything that needs to be there in orderfor the modules of the chip to perform their functions into a memory onto the chip or attached to the chip.

Additionally, the API may be configured to allow the chip to interfacewith the circuit board of the sequencer, when included therewith, so asto receive the FASTQ sequencing files directly from the sequencer suchas immediately once they have been generated and then transfers thatinformation to the configuration manager which then directs thatinformation to the appropriate memory banks in the hardware that makesthat information available to the pertinent modules of the hardware sothat they can perform their designated functions on that information soas to call bases, map, align, sort, etc. the sample DNA with respect tothe reference genome.

Further still, a client level interface (CLI) may be included whereinthe CLI may allow the user to call one or more of these functionsdirectly. In various embodiments, the CLI may be a software applicationthat is adapted to configure the use of the hardware. The CLI,therefore, may be a program that accepts instructions, e.g., arguments,and makes functionality available simply by calling an applicationprogram. As indicated above, the CLI can be command line based or GUI(graphical user interface) based. The line based commands happen at alevel below the GUI, where the GUI includes a windows based file managerwith click on function boxes that delineate which modules will be usedand the parameters of their use. For example, in operation, ifinstructed, the CLI will locate the reference, will determine if a hashtable and/or index needs to be generated, or if already generated locatewhere it is stored, and direct the uploading of the generated hash tableand/or index, etc. These type of instructions may appear as user optionsat the GUI that the user can select the chip to perform.

Furthermore, a library may be included wherein the library may includepre-existing, editable, configuration files, such as files orientated tothe typical user selected functioning of the hardware, such as withrespect to a portion or whole genome analysis, for instance, forancestry analysis, or disease diagnostics, or drug discovery, or proteinprofiling, etc. These types of preset parameters, such as for performingsuch analyses, may be stored in the library. For example, if theplatform herein described is employed such as for oncology research, thepreset parameters may be configured differently than if the platformwere directed simply to researching a genealogy.

More particularly, for oncology, accuracy may be an important factor,therefore, the parameters of the system may be set to ensure increasedaccuracy albeit in exchange for possibly a decrease in speed. However,for other genomics applications, speed may be the key determinant andtherefore the parameters of the system may be set to maximize speed,which however may sacrifice some accuracy. Accordingly, in variousembodiments, often used parameter settings for performing differenttasks can be preset into the library to facilitate ease of use. Suchparameter settings may also include the necessary software applicationsemployed in running the system. For instance, the library may containthe code that executes the API, and may further include sample files,scripts, and any other ancillary information necessary for running thesystem. Hence, the library may be configured for compiling software forrunning the API as well as various executables.

In various instances, the chip may also include a memory, such as aRandom Access Memory (RAM) or a Dynamic Rapid Access Memory with e.g. aDDR3 interface, such as a memory that may be used for facilitating theperformance of the various modules described herein, for instance, themapper, aligner, and/or sorter. For example, the DRAM may be where thereference, the hash table, and/or the hash table index, and/or reads maybe stored. Further, the memory may be used for facilitating theperformance of various other modules described herein, for instance, thededuper, local realigner, base quality score recalibrator, variantcaller, compressor, and/or decompresor. For example, the DRAM may bewhere sorted reads, annotated reads, compressed reads, and/or variantcalls may be stored. Further, the memory may be configured so as toinclude a separate interface for each of the various memory modulesemployed by the aligner and/or any other module, such as where eachmemory may include a file layer and logical layer. As indicated above,because there may be multiple memories and/or multiple modules, a chiplevel interconnect may be included so as to facilitate communicationthrough the chip.

Accordingly, in various instances, an apparatus of the disclosure mayinclude a chip, wherein the chip includes an integrated circuit that isformed of a set of hardwired digital logic circuits that may beinterconnected by one or more physical electrical interconnects. Invarious embodiments, the one or more physical electrical interconnectsinclude an input to the integrated circuit that may be connected with anelectronic data source for receiving data. Further, in certainembodiments, the hardwired digital logic circuits may be arranged as aset of processing engines, such as wherein each processing engine may beformed of a subset of the hardwired digital logic circuits, which areconfigured to perform one or more of the steps in the sequence analysispipeline. More particularly, each subset of the hardwired digital logiccircuits may be in a wired configuration so as to perform the one ormore steps in the sequence analysis pipeline.

In various instances, the set of processing engines may include one ormore of a mapping module, an alignment module, and/or a sorting module,such as where the one or more of these modules are in the wiredconfiguration. For instance, a mapping module may be included, where inthe wired configuration, the mapping module may access an index, such asof one or more genetic reference sequences, e.g., from a memory, such asvia one or more of the plurality of physical electronic interconnects,so as to map the plurality of reads to one or more segments of the oneor more genetic reference sequences. Further, in various instances, analignment module may be included, wherein the wired configuration, thealignment module may access the one or more genetic reference sequences,e.g., from the memory, such as via one or more of the plurality ofphysical electronic interconnects, so to align the plurality of reads tothe one or more segments of the one or more genetic reference sequences.Further still, in various instances, a sorting module may be included,wherein the wired configuration, the sorting module may access the oneor more aligned sequences, e.g., from the memory, such as via one ormore of the plurality of physical electronic interconnects, so to sortthe plurality of reads to a chromosome, such as from the one or moregenetic reference sequences. In like manner, in various instances, oneor more of local realignment, duplicate marking, base quality scorerecalibration, and/or variant calling modules may be included in thechip, such as in the wired configuration consistent as with the modulesdescribed above, so as to perform their respective functions.

Further, as indicated above, in various instances a chip of thedisclosure may be configured as an expansion card, such as where thechip includes a PCIe bus and is positioned so as to be in communicationwith one or more memories, such as being surrounding by memories, suchas being substantially surrounded by memories, such as being entirelysurrounded by memories. In various embodiments, the chip may be a denseand/or fast FPGA chip, that in various instances, may be convertible toan ASIC. As indicated above, the modules herein disclosed may beimplemented in the hardware of the chip, such as by being hardwiredtherein, and in such instances their implementation may be such thattheir functioning may take place at a faster speed as compared to whenimplemented in software, such as where there are minimal instructions tobe fetched, read, and/or executed. Hence, given the unique hardwareimplementation, the modules of the disclosure may function directly inaccordance with their operations parameters, such as without needing tofetch, read, and/or execute instructions. Additionally, memoryrequirements and processing times may be reduced, such as where thecommunications within chip is via files rather than through accessing amemory. Of course, in some instances, the chip and/or card may be sizedso as to include more memory, such as more on board memory, so as toenhance parallel processing capabilities, thereby resulting in evenfaster processing speeds. For instance, in certain embodiments, a chipof the disclosure may include an embedded DRAM, so that the chip doesnot have to rely on external memory, which would therefore result in afurther increase in processing speed, such as where a Burrows-Wheeleralgorithm may be employed, instead of a hash table and hash function,which may in various instances, rely on external, e.g., host memory. Insuch instances, the running of the entire pipeline can be accomplishedin 6 minutes or less, such as from start to finish.

As indicated above, there are various different points where any givenmodule can be positioned on the hardware, or be positioned remotelytherefrom, such as on a server accessible on the cloud. Where a givenmodule is positioned on the chip, e.g., hardwired into the chip, itsfunction may be performed by the hardware, however, where desired, themodule may be positioned remotely from the chip, at which point theplatform may include the necessary instrumentality for sending therelevant data to a remote location, such as a server accessible via thecloud, so that the particular module's functionality may be engaged forfurther processing of the data, in accordance with the user selecteddesired protocols. Accordingly, part of the platform may include aweb-based interface for the performance of one or more tasks pursuant tothe functioning of one or more of the modules disclosed herein. Forinstance, where mapping, alignment, and/or sorting are all modules thatmay occur on the chip, in various instances, one or more of localrealignment, duplicate marking, base quality core recalibration, and/orvariant calling may take place on the cloud.

Additionally, in various embodiments, all of mapping, aligning, andsorting, may take place on the chip, and local realignment, duplicatemarking, and/or base quality score recalibration may, in variousembodiments, also take place on the chip, and in various instances,various compression protocols, such as BAM and CRAM, may also take placeon the chip. However, once the data is compressed it may be sent up tothe cloud, such as for the performance of the variant calling module.This might be useful especially given the fact that variant calling canbe a moving target, e.g., there is not one standardized agreed uponalgorithm that the industry uses. Hence, different algorithms can beemployed to achieve a different type of result, and as such having acloud based module for the performance of this function may be usefulfor allowing the flexibility to select which algorithm is useful at anyparticular given moment, and also as for serial and/or parallelprocessing. Accordingly, any one of the modules disclosed herein can beimplemented as either hardware, e.g., on the chip, or software, e.g., onthe cloud, but in certain embodiments, all of the modules may beconfigured so that their function may be performed on the chip, or allof the modules may be configured so that their function may be performedremotely, such as on the cloud, or there will be a mixture of moduleswherein some are positioned on the chip and some are positioned on thecloud. Further, as indicated, in various embodiments, the chip itselfmay be configured so as to function in conjunction with, and in someembodiments, in immediate operation with a genetic sequencer.

More specifically, in various embodiments, an apparatus of thedisclosure may be a chip, such as a chip that is configured forprocessing genomics data, such as by employing a pipeline of dataanalysis modules. According, as can be seen with respect to FIG. 1, agenomics pipeline processor chip 100 is provided along with associatedhardware of a genomics pipeline processor system 10. The chip 100 hasone or more connections to external memory 102 (at “DDR3 MemController”), and a connection 104 (e.g., “PCIe Interface”) to theoutside world, such as a host computer 106, for example. A crossbar 108(e.g., switch) provides access to the memory interfaces to variousrequestors. DMA engines 110 transfer data at high speeds between thehost and the processor chip's 100 external memories 102 (via thecrossbar 108), and/or between the host and a central controller 112. Thecentral controller 112 controls chip operations, especially coordinatingthe efforts of multiple processing engines. The processing engines areformed of a set of hardwired digital logic circuits that areinterconnected by physical electrical interconnects, and are organizedinto engine clusters 114. In some implementations, the engines in onecluster share one crossbar port, via an arbiter. The central controller112 has connections to each of the engine clusters. Each engine cluster114 has a number of processing engines for processing genomic data,including a mapper 120 (or mapping module), an aligner 122 (or aligningmodule), and a sorter 124 (or sorting module). An engine cluster 114 caninclude other engines or modules, as well.

In accordance with one data flow model consistent with implementationsdescribed herein, the host sends commands and data via the DMA engines110 to the central controller 112, which load-balances the data to theprocessing engines. The processing engines return processed data to thecentral controller 112, which streams it back to the host via the DMAengines 110. This data flow model is suited for mapping and alignment.

In accordance with an alternative data flow model consistent withimplementations described herein, the host streams data into theexternal memory, either directly via DMA engines 110 and the crossbar108, or via the central controller 112. The host sends commands to thecentral controller 112, which sends commands to the processing engines,which instruct the processing engines as to what data to process. Theprocessing engines access input data from the external memory, processit, and write results back to the external memory, reporting status tothe central controller 112. The central controller 112 either streamsthe result data back to the host from the external memory, or notifiesthe host to fetch the result data itself via the DMA engines 110.

FIG. 2 illustrates a genomics pipeline processor system 20, showing afull complement of processing engines inside an engine cluster 214. Thepipeline processor system 20 may include one or more engine clusters214. In some implementations, the pipeline processor system 20 includesfour our more engine clusters 214. The processing engines or processingengine types can include, without limitation, a mapper, an aligner, asorter, a local realigner, a base quality recalibrater, a duplicatemarker, a variant caller, a compressor and/or a decompressor. In someimplementations, each engine cluster 214 has one of each processingengine type. Accordingly, all processing engines of the same type canaccess the crossbar 208 simultaneously, through different crossbarports, because they are each in a different engine cluster 214. Notevery processing engine type needs to be formed in every engine cluster214. Processing engine types that require massive parallel processing ormemory bandwidth, such as the mapper (and attached aligner(s)) andsorter, may appear in every engine cluster of the pipeline processorsystem 20. Other engine types may appear in only one or some of theengine clusters 214, as needed to satisfy their performance requirementsor the performance requirements of the pipeline processor system 20.

FIG. 3 illustrates a genomics pipeline processor system 30, showing, inaddition to the engine clusters described above, one or more embeddedcentral processing units (CPUs) 302. Examples of such embedded CPUsinclude Snapdragons® or standard ARM® cores. These CPUs execute fullyprogrammable bio-IT algorithms, such as advanced variant calling. Suchprocessing is accelerated by computing functions in the engine clusters,which can be called by the CPU cores 302 as needed. Furthermore, evenengine-centric processing, such as mapping and alignment, can be managedby the CPU cores 302, giving them heightened programmability.

FIG. 4 illustrates a processing flow for a genomics pipeline processorsystem and method. In some preferred implementations, there are threepasses over the data. The first pass includes mapping 402 and alignment404, with the full set of reads streamed through the engines. The secondpass includes sorting 406, where one large block to be sorted (e.g., asubstantial portion or all reads previously mapped to a singlechromosome) is loaded into memory, sorted by the processing engines, andreturned to the host. The third pass includes downstream stages (localrealignment 408, duplicate marking 410, base quality score recalibration(BQSR) 412, BAM output 414, reduced BAM output 416, and/or CRAMcompression 418). The steps and functions of the third pass may be donein any combination or subcombination, and in any order, in a singlepass. A virtual pipeline architecture, such as described above, is usedto stream reads from the host into circular buffers in memory, throughone processing engine after another in sequence, and back out to thehost. In some implementations, CRAM decompression can be a separatestreaming function. In some implementations, the BAM output 414, reducedBAM output 416, and/or CRAM compression 418 can be replaced with variantcalling, compression and decompression.

FIG. 5 shows a general block diagram of the current invention. In Block1 a hardware implementation of a sequence analysis pipeline isdescribed. This can be done in a number of different ways such as anFPGA or ASIC implementation. The functional blocks that are implementedby the FPGA or ASIC are shown in FIG. 5. FIG. 5 includes a number ofblocks or modules to do sequence analysis. The input to the hardwarerealization can be a FASTQ file, but is not limited to this format. Inaddition to the FASTQ file, the input to the FPGA or ASIC consists ofside information, such as Flow Space Information from technology such asthe Ion Torrent. The blocks or modules in FIG. 5 illustrate thefollowing blocks: Error Control, Mapping, Alignment, Sorting, LocalRealignment, Duplicate Marking, Base Quality Recalibration, BAM and SideInformation reduction and variant calling.

These blocks or modules can be present inside, or implemented by, thehardware, but some of these blocks may be omitted or other blocks addedto achieve the purpose of realizing a sequence analysis pipeline. Blocks2 and 3 describe two alternatives of a The sequence analysis pipelineplatform. The sequence analysis pipeline platform comprising an FPGA orASIC and software assisted by a host (i.e., PC, server, cluster or cloudcomputing) with cloud and/or cluster storage. Blocks 4-7 describedifferent interfaces that the sequence analysis pipeline can have. InBlocks 4 and 6 the interface can be a PCIe interface, but is not limitedto a PCIe interface. In Blocks 5 and 7 the hardware (FPGA or ASIC) canbe directly integrated into a sequencing machine. Blocks 8 and 9describe the integration of the hardware sequence analysis pipelineintegrated into a host system such as a PC, server cluster or sequencer.Surrounding the hardware FPGA or ASIC are lots of DDR3 memory elementsand a PCIe interface. The board with the FPGA/ASIC connects to a hostcomputer, consisting of a host CPU, that could be either a low power CPUsuch as an ARM®, Snapdragon®, or any other processor. Block 10illustrates a hardware sequence analysis pipeline API that can beaccessed by third party applications to perform tertiary analysis.

Accordingly, in various embodiments, an apparatus of the disclosure mayinclude a computing architecture, such as embedded in a siliconapplication specific integrated circuit (ASIC) 100 as seen in FIGS. 6and 7. The ASIC 100 can be integrated into a printed circuit board (PCB)104, such as a Peripheral Component Interface—Express (PCIe) card, thatcan be plugged into a computing platform. In various instances, as shownin FIG. 6, the PCIe card 104 may include a single ASIC 100, which ASICmay be surrounded by local memories 105, however, in variousembodiments, the PCIe card 104 may include a plurality of ASICs 100A,100B and 100C. In various instances, the PCI card may also include aPCIe bus. This PCIe card 104 can be added to a computing platform toexecute algorithms on extremely large data sets. Accordingly, in variousinstances, the overall work flow of genomic sequencing involving theASIC may include the following: Sample preparation, Alignment (includingmapping and alignment), Variant analysis, Biological Interpretation,and/or Specific Applications.

Hence, in various embodiments, an apparatus of the disclosure mayinclude a computing architecture that achieves the high performanceexecution of algorithms, such as mapping and alignment algorithms, thatoperate on extremely large data sets, such as where the data setsexhibit poor locality of reference (LOR). These algorithms are designedto reconstruct a whole genome from millions of short read sequences,from modern so-called next generation sequencers, require multi-gigabytedata structures that are randomly accessed. Once reconstruction isachieved, as described herein above, further algorithms with similarcharacteristics are used to compare one genome to libraries of others,do gene function analysis, etc.

Currently, there are two major approaches in use, general purposemulticore CPUs and general purpose Graphic Processing Units (GPGPUs). Insuch an instance ach CPU in a multicore system may have a classicalcache based architecture, wherein instructions and data are fetched froma level 1 cache (L1 cache) that is small but has extremely fast access.Multiple L1 caches may be connected to a larger but slower shared L2cache. The L2 cache may be connected to a large but slower DRAM (DynamicRandom Access Memory) system memory, or may be connected to an evenlarger but slower L3 cache which may then connected to DRAM. Anadvantage of this arrangement may be that applications in which programsand data exhibit locality of reference behave nearly as if they areexecuting on a computer with a single memory as large as the DRAM but asfast as the L1 cache. Because full custom, highly optimized CPUs operateat very high clock rates, e.g., 2 to 4 GHz, this architecture may beessential to achieving good performance.

Further, GPGPUs may be employed to extend this architecture, such as byimplementing very large numbers of small CPUs, each with their own smallL1 cache, wherein each CPU executes the same instructions on differentsubsets of the data. This is a so called SIMD (Single Instructionstream, Multiple Data stream) architecture. Economy is gained by sharingthe instruction fetch and decode logic across a large number of CPUs.Each cache has access to multiple large external DRAMs via aninterconnection network. Assuming the computation to be performed ishighly parallelizable, GPGPUs have a significant advantage over generalpurpose CPUs due to having large numbers of computing resources.Nevertheless, they still have a caching architecture and theirperformance is hurt by applications that do not have a high enoughdegree of locality of reference. That leads to a high cache miss rateand processors that are idle while waiting for data to arrive from theexternal DRAM.

For instance, in various instances, Dynamic RAMs may be used for systemmemory because they are more economical than Static RAMs (SRAM). Therule of thumb used to be that DRAMs had 4× the capacity for the samecost as SRAMs. However, due to declining demand for SRAMs in favor ofDRAMs, that difference has increased considerably due to the economiesof scale that favor DRAMs which are in high demand. Independent of cost,DRAMs are 4× as dense as SRAMs laid out in the same silicon area becausethey only require one transistor and capacitor per bit compared to 4transistors per bit to implement the SRAM's flip-flop. The DRAMrepresents a single bit of information as the presence or absence ofcharge on a capacitor. A problem with this arrangement is that thecharge decays over time, so it has to be refreshed periodically. Theneed to do this has led to architectures that organize the memory intoindependent blocks and access mechanisms that deliver multiple words ofmemory per request. This compensates for times when a given block isunavailable while being refreshed. The idea is to move a lot of datawhile a given block is available. This is in contrast to SRAMs in whichany location in memory is available in a single access in a constantamount of time. This characteristic allows memory accesses to be singleword oriented rather than block oriented. DRAMs work well in a cachingarchitecture because each cache miss leads to a block of memory beingread in from the DRAM. The theory of locality of reference is that ifjust accessed word N, then probably going to access words N+1, N+2, N+3and so on, soon.

FIG. 8 illustrates a system 500 for executing a sequence analysispipeline on genetic sequence data. The system 500 includes aconfiguration manager 502 that includes a computing system. Thecomputing system of the configuration manager 502 can include a personalcomputer or other computer workstation, or can be implemented by a suiteof networked computers. The configuration manager 502 can furtherinclude one or more third party applications connected with thecomputing system by one or more APIs, which, with one or moreproprietary applications, generate a configuration for processinggenomics data from a sequencer or other genomics data source. Theconfiguration manager 502 further includes drivers that load theconfiguration to the genomics pipeline processor system 10. The genomicspipeline processor system 10 can output result data to, or be accessedvia, the Web 504 or other network, for storage of the result data in anelectronic health record 506 or other knowledge database 508.

In some implementations, the chip implementing the genomics pipelineprocessor can be connected or integrated in a sequencer. The chip canalso be connected or integrated on an expansion card, e.g. PCIe, and theexpansion card can by connected or integrated in a sequencer. In otherimplementations, the chip can be connected or integrated in a servercomputer that is connected to a sequencer, to transfer genomic readsfrom the sequencer to the server. In yet other implementations, the chipcan be connected or integrated in a server in a cloud computing clusterof computers and servers. A system can include one or more sequencersconnected (e.g. via Ethernet) to a server containing the chip, wheregenomic reads are generated by the multiple sequencers, transmitted tothe server, and then mapped and aligned in the chip.

The memory architecture can consist of M memory modules that interfacewith an ASIC. The ASIC may be implemented using many differenttechnologies, including FPGAs (Field Programmable Gate Arrays), standardcells, or full custom logic. Within the ASIC are a Memory Subsystem(MSS) and Functional Processing Units (FPUs). The MSS contains M memorycontrollers (MCs) for the memory modules, N system memory interfaces(SMIs) for the FPUs, and an N×M crossbar that allows any SMI to accessany MC. Arbitration is provided in the case of contention.

Each memory module is constructed from DRAM chips that are addressed byan AMM bit word and support data transfers D_(MM) bits wide. The memoryhas 2^(A) ^(MM) address locations. A key characteristic of DRAM is thatit performs reads/writes in W word bursts using the supplied address asthe base address, B, and fetching or storing locations B+1, B+2, . . .B+W−1 as well. A typical value for W is 8.

In the MSS of the ASIC, each memory controller supplies the requiredcontrol signals and performs any necessary multiplexing/demultiplexingbetween the system word width, D_(SYS), and the memory word width,D_(MM), as well as handling the requirements for read/write bursts. Itcan contain extra buffering so that multiple memory requests can bequeued up and processed in a pipelined fashion to maximize throughput.This compensates for multiple clock cycles of latency betweenpresentation of an address and completion of a memory operation (read orwrite).

The MC necessarily operates at the speed of the attached DRAM in amemory module. Assume its clock rate is C_(MM). This is often severaltimes faster than the core speed at which the majority of the logic inthe ASIC operates which is C_(SYS). Hence themultiplexing/demultiplexing logic is placed close to its associatedinterface pins to minimize signal distances. Demultiplexing is the firstoperation performed on incoming data and multiplexing is the lastoperation performed on outgoing data. The remainder of the MSS operateson D_(SYS) width data which is wider than D_(MM), enabling use of theslower C_(SYS) clock speed.

Each system memory interface in the MSS presents an A_(SYS) bit addressbus and a D_(SYS) bit data bus to any attached FPU. The SMI is designedto make it appear to an attached FPU that it has random access to asingle large fast memory. The FPU has no awareness of the existence ofseparate memory modules. A_(SYS) is large enough to allow access to anymemory location in any attached memory module. The mapping from systemaddress space to memory module address space is explained below.

The N system memory interfaces are cross connected to the M memorymodules via an N×M crossbar. The crossbar provides min(M,N) simultaneousconnections among the SMIs and MCs, provides arbitration for conflicts,and facilitates translation of system address space into memory moduleaddress space.

The organization of FPUs is highly flexible. One or more FPUs can sharethe same system memory interface. To maximize performance, FPUs that donot operate at the same time should share an SMI. Those that operateconcurrently, should be attached to different SMIs. An FPU that operateson a data structure larger than D_(SYS) can use multiple SMIs to accessthe whole data structure in a single memory operation. Hence this memoryarchitecture supports a wide range of computation architectures. EachFPU may be identical and thus an array of them may be implemented in atwo dimensional structure. This is illustrated in Error! Referencesource not found. where FPU(i,j) is the j^(th) unit attached to SMI i,0<i<N, 0<j<k_(i). In this case, all the k_(i) are the same size andk_(i) may be as small as 1. This supports SIMD (single instructionstream, multiple data stream) and MIMD architectures (multipleinstruction stream, multiple data stream) depending on whether the FPUsfetch instructions from the same or individual instruction memories.

One or more aspects or features of the subject matter described hereincan be realized in digital electronic circuitry, integrated circuitry,specially designed application specific integrated circuits (ASICs),field programmable gate arrays (FPGAs) computer hardware, firmware,software, and/or combinations thereof.

These various aspects or features can include implementation in one ormore computer programs that are executable and/or interpretable on aprogrammable system including at least one programmable processor, whichcan be special or general purpose, coupled to receive data andinstructions from, and to transmit data and instructions to, a storagesystem, at least one input device, and at least one output device. Theprogrammable system or computing system may include clients and servers.A client and server are generally remote from each other and typicallyinteract through a communication network. The relationship of client andserver arises by virtue of computer programs running on the respectivecomputers and having a client-server relationship to each other.

These computer programs, which can also be referred to as programs,software, software applications, applications, components, or code,include machine instructions for a programmable processor, and can beimplemented in a high-level procedural and/or object-orientedprogramming language, and/or in assembly/machine language. As usedherein, the term “machine-readable medium” refers to any computerprogram product, apparatus and/or device, such as for example magneticdiscs, optical disks, memory, and Programmable Logic Devices (PLDs),used to provide machine instructions and/or data to a programmableprocessor, including a machine-readable medium that receives machineinstructions as a machine-readable signal. The term “machine-readablesignal” refers to any signal used to provide machine instructions and/ordata to a programmable processor. The machine-readable medium can storesuch machine instructions non-transitorily, such as for example as woulda non-transient solid-state memory or a magnetic hard drive or anyequivalent storage medium. The machine-readable medium can alternativelyor additionally store such machine instructions in a transient manner,such as for example as would a processor cache or other random accessmemory associated with one or more physical processor cores.

To provide for interaction with a user, one or more aspects or featuresof the subject matter described herein can be implemented on a computerhaving a display device, such as for example a cathode ray tube (CRT), aliquid crystal display (LCD) or a light emitting diode (LED) monitor fordisplaying information to the user and a keyboard and a pointing device,such as for example a mouse or a trackball, by which the user mayprovide input to the computer. Other kinds of devices can be used toprovide for interaction with a user as well. For example, feedbackprovided to the user can be any form of sensory feedback, such as forexample visual feedback, auditory feedback, or tactile feedback; andinput from the user may be received in any form, including, but notlimited to, acoustic, speech, or tactile input. Other possible inputdevices include, but are not limited to, touch screens or othertouch-sensitive devices such as single or multi-point resistive orcapacitive trackpads, voice recognition hardware and software, opticalscanners, optical pointers, digital image capture devices and associatedinterpretation software, and the like.

The subject matter described herein can be embodied in systems,apparatus, methods, and/or articles depending on the desiredconfiguration. The implementations set forth in the foregoingdescription do not represent all implementations consistent with thesubject matter described herein. Instead, they are merely some examplesconsistent with aspects related to the described subject matter.Although a few variations have been described in detail above, othermodifications or additions are possible. In particular, further featuresand/or variations can be provided in addition to those set forth herein.For example, the implementations described above can be directed tovarious combinations and subcombinations of the disclosed featuresand/or combinations and subcombinations of several further featuresdisclosed above. In addition, the logic flows depicted in theaccompanying figures and/or described herein do not necessarily requirethe particular order shown, or sequential order, to achieve desirableresults. Other implementations may be within the scope of the followingclaims.

What is claimed is:
 1. A system for executing a sequence analysispipeline on genetic sequence data, the system comprising: a memory forstoring one or more genetic reference sequences, an index of the one ormore genetic reference sequences, and a plurality of reads of genomicdata, each of the genetic reference sequences and the plurality of readsof genomic data comprising a sequence of nucleotides; an integratedcircuit formed of a set of hardwired digital logic circuits that areinterconnected by a plurality of physical electrical interconnects, oneor more of the plurality of physical electrical interconnects comprisinga memory interface for the integrated circuit to access the memory, thehardwired digital logic circuits being arranged as a set of processingengines to perform steps in the sequence analysis pipeline on theplurality of reads of genomic data, the set of processing enginescomprising: a mapping module in a first configuration of the hardwireddigital logic circuits, the mapping module comprising: a first subset ofthe first configuration of the hardwired digital logic circuits that isconfigured to receive a read of genomic data via one or more of theplurality of physical electrical interconnects; a second subset of thefirst configuration of the hardwired digital logic circuits that isconfigured to extract a portion of the read to generate a seed, the seedrepresenting a subset of the sequence of nucleotides represented by theread; a third subset of the first configuration of the hardwired digitallogic circuits that is configured to calculate an address within theindex based on the seed; a fourth subset of the first configuration ofthe hardwired digital logic circuits that is configured to access theaddress in the index in the memory; a fifth subset of the firstconfiguration of the hardwired digital logic circuits that is configuredto receive a record from the address, the record representing positioninformation in the genetic reference sequence; and a sixth subset of thefirst configuration of the hardwired digital logic circuits that isconfigured to determine one or more matching positions from the read tothe genetic reference sequence based on the record, so as to producemapped positions for the read; and an alignment module in a secondconfiguration of the hardwired digital logic circuits, the alignmentmodule comprising: a first subset of the second configuration of thehardwired digital logic circuits that is configured to receive one ormore of the mapped positions for the read from the memory; a secondsubset of the second configuration of the hardwired digital logiccircuits that is configured to access the memory to retrieve a segmentof the genetic reference sequence corresponding to each of the mappedpositions; a third subset of the second configuration of the hardwireddigital logic circuits that is configured to calculate an alignment ofthe read to each retrieved segment of the genetic reference sequence; afourth subset of the second configuration of the hardwired digital logiccircuits that is configured to generate a score for each alignment; anda fifth subset of the second configuration of the hardwired digitallogic circuits that is configured to select at least one best-scoringalignment of the read; and an output formed of a subset of the pluralityof physical electrical interconnects for outputting the best-scoringalignment of the read.
 2. The system according to claim 1, wherein theintegrated circuit is a field programmable gate array (FPGA).
 3. Thesystem according to claim 2, wherein the set of processing engines ofthe integrated circuit further comprises a sorting module in anotherhardwired configuration to sort each aligned read according to the oneor more positions in the genetic reference sequence.
 4. The systemaccording to claim 2, wherein the set of processing engines of theintegrated circuit further comprises a variant call module in anotherhardwired configuration to perform a variant call operation on eachaligned read according to the one or more positions in the geneticreference sequence.
 5. The system according to claim 1, wherein thealignment module is configured for performing a gapless alignment. 6.The system according to claim 1, wherein the alignment module is furtherconfigured for performing a Smith Waterman alignment.
 7. A system forexecuting a sequence analysis pipeline on genetic sequence data, thesystem comprising: a memory for storing one or more genetic referencesequences, an index of the one or more genetic reference sequences, anda plurality of reads of genomic data, each of the genetic referencesequences and the plurality of reads of genomic data comprising asequence of nucleotides; and an integrated circuit formed of a set ofhardwired digital logic circuits that are interconnected by a pluralityof physical electrical interconnects, one or more of the plurality ofphysical electrical interconnects comprising a memory interface for theintegrated circuit to access the memory, the hardwired digital logiccircuits being arranged as a set of processing engines to perform stepsin the sequence analysis pipeline on the plurality of reads of genomicdata, the set of processing engines comprising: a mapping module in afirst configuration of the hardwired digital logic circuits, the mappingmodule comprising: a first subset of the first configuration of thehardwired digital logic circuits that is configured to receive a read ofgenomic data via one or more of the plurality of physical electricalinterconnects; a second subset of the first configuration of thehardwired digital logic circuits that is configured to extract a portionof the read to generate a seed, the seed representing a subset of thesequence of nucleotides represented by the read; a third subset of thefirst configuration of the hardwired digital logic circuits that isconfigured to calculate an address within the index based on the seed; afourth subset of the first configuration of the hardwired digital logiccircuits that is configured to access the address in the index in thememory; a fifth subset of the first configuration of the hardwireddigital logic circuits that is configured to receive a record from theaddress, the record representing position information in the geneticreference sequence; and a sixth subset of the first configuration of thehardwired digital logic circuits that is configured to determine one ormore matching positions from the read to the genetic reference sequencebased on the record; and an alignment module in a second configurationof the hardwired digital logic circuits to access the one or moregenetic reference sequences from the memory via the memory interface toalign the mapped read from the mapping module to one or more positionsin the one or more segments of the one or more genetic referencesequences to produce aligned read data; and an output formed of a subsetof the plurality of physical electrical interconnects for outputting thealigned read data.
 8. The system according to claim 7, wherein theintegrated circuit is a field programmable gate array (FPGA).
 9. Thesystem according to claim 8, wherein the set of processing engines ofthe integrated circuit further comprises a sorting module in anotherhardwired configuration to sort each aligned read according to the oneor more positions in the genetic reference sequence.
 10. The systemaccording to claim 8, wherein the set of processing engines of theintegrated circuit further comprises a variant call module in anotherhardwired configuration to perform a variant call operation on eachaligned read according to the one or more positions in the geneticreference sequence.
 11. The system according to claim 7, wherein thealignment module is configured for performing a gapless alignment. 12.The system according to claim 7, wherein the alignment module is furtherconfigured for performing a Smith Waterman alignment.
 13. The systemaccording to claim 7, wherein the index of genetic reference datafurther comprises a hash table, and wherein the mapping module applies ahash function to the at least some of the sequence of nucleotides toaccess the hash table of the index.
 14. A system for executing asequence analysis pipeline on genetic sequence data, the systemcomprising: a memory for storing one or more genetic referencesequences, an index of the one or more genetic reference sequences, anda plurality of reads of genomic data, each of the genetic referencesequences and the plurality of reads of genomic data comprising asequence of nucleotides; and an integrated circuit formed of a set ofhardwired digital logic circuits that are interconnected by a pluralityof physical electrical interconnects, one or more of the plurality ofphysical electrical interconnects comprising a memory interface for theintegrated circuit to access the memory, the hardwired digital logiccircuits being arranged as a set of processing engines to perform stepsin the sequence analysis pipeline on the plurality of reads of genomicdata, the set of processing engines comprising: a mapping module in afirst configuration of the hardwired digital logic circuits to access inthe memory, via the memory interface, at least some of the sequence ofnucleotides in a selected read of the plurality of reads and the indexof the one or more genetic reference sequences, and to map the selectedread to one or more segments of the one or more genetic referencesequences based on the index to produce a mapped read; and an alignmentmodule in a second configuration of the hardwired digital logiccircuits, the alignment module comprising: a first subset of the secondconfiguration of the hardwired digital logic circuits that is configuredto receive one or more of mapped reads from the memory; a second subsetof the second configuration of the hardwired digital logic circuits thatis configured to access the memory to retrieve a segment of the geneticreference sequence corresponding to each of the received mapped reads; athird subset of the second configuration of the hardwired digital logiccircuits that is configured to calculate an alignment of the read toeach retrieved segment of the genetic reference sequence; a fourthsubset of the second configuration of the hardwired digital logiccircuits that is configured to generate a score for each alignment; anda fifth subset of the second configuration of the hardwired digitallogic circuits that is configured to select at least one best-scoringalignment of the read; and an output formed of a subset of the pluralityof physical electrical interconnects for outputting the best-scoringalignment of the read.
 15. The system according to claim 14, wherein theintegrated circuit is a field programmable gate array (FPGA).
 16. Thesystem according to claim 15, wherein the set of processing engines ofthe integrated circuit further comprises a sorting module in anotherhardwired configuration to sort each aligned read according to the oneor more positions in the genetic reference sequence.
 17. The systemaccording to claim 15, wherein the set of processing engines of theintegrated circuit further comprises a variant call module in anotherhardwired configuration to perform a variant call operation on eachaligned read according to the one or more positions in the geneticreference sequence.
 18. The system according to claim 14, wherein thealignment module is configured for performing a gapless alignment. 19.The system according to claim 14, wherein the alignment module isfurther configured for performing a Smith Waterman alignment.
 20. Thesystem according to claim 14, wherein the index of genetic referencedata further comprises a hash table, and wherein the mapping moduleapplies a hash function to the at least some of the sequence ofnucleotides to access the hash table of the index.